I am trying to mathematically model the temperature distribution and evolution in a greenhouse, which conserves heat due to the greenhouse effect. Here is a transverse schematic ($H$ is the height of the box, $l$ its width and $L$ is length but is not represented in this 2D schema) of the situation I'm trying to solve:
There are several radiative energy flows (visible $F_d$ from the sun, infrared $F_i$ from the atmosphere, infrared $F_s=\sigma T_s^4$ from the black soil and infrared $F_p=\sigma T_p^4$ from the plastic/glass wall). The plastic/glass absorbs infrared but not visible light, it is why $F_d$ goes directly into the black soil in the schema. There are also convective flows inside the greenhouse due to the airflow $Q$ but also due to natural convection between the two plates of different temperature ($T_s > T_p$). Adimensionnal numbers and correlations can be used to deduce the convective factor $h$.
For the moment, I made the hypothesis that $F_d$ and $F_i$ were constant. I also have considered a stationary situation where the energy does not accumulate, the outgoing energy is thus equal to the incoming energy. By doing such energy balance locally (air and soil) and globally, I have obtained the following algebraic system of non-linear equation, with three unknows $P$, $T_s$ and $T_p$, where $P$ is the heat given to the air inside the greenhouse:
\begin{array}{r c l} P = h(T_p-T) + h(T_s-T) \\ F_d + F_i = P + \sigma T_p^4 \\ F_d + \sigma T_p^4 = \sigma T_s^4 + h(T_s-T) \end{array}
Using the Newton-Raphson method, I was able to solve this non-linear system, but I would now like to refine my mathematical model. So here is my two problems:
$F_d(t)=A\sin(\pi (t-t_r)/(t_s-t_r))$ is not constant (where $t_r$ and $t_s$ are respectively the sunrise and sunset time). The hypothesis of a stationary problem is therefore erroneous. I tried to modify my equation system but I'm stuck because I don't know where to put the time derivatives and the appropriate constants (heat capacity, volumic mass ?). Once I have the differential equation system, I think I know how to proceed, using Crank-Nicolson discretization method and returning to an algebraic system that can be solved by Newton-Raphson.
In addition to being able to have the temporal evolution of temperatures $T$ (air), $T_s$ (black soil) and $T_p$ (plastic/glass plate), I would also like to have their spatial distribution. What I mean is the value of $T(z,t)$ for all position $z$ and for all time $t$, knowing the initial condition $T(z,t_0)$ for all $z$, and assuming a constant air flow $Q$. I guess because of the airflow $Q$, there will be a air transport equation involved, but I'm not sure how to include these new equations.
I will be infinitely grateful to anyone who can help. It is a personal project that I have set myself for educational purposes, I am not an expert but I try as best I can to carry out this project.