#. # Parallelepiped whose Young’s modulus is a function of the temperature
#.
#. The thermo-mechanical problem is solved in two stages. First, the heat conduction
#. equation is solved over a coarse first-order mesh to find the non-dimensional temperature distribution.
#. Then, a mechanical problem is solved where $T(x,y,z)$ is read from the first mesh
#. and interpolated into a finer second-order mesh so to as evaluate the non-dimensional Young's modulus
#. as
#.
#. $$E\Big(T(x,y,z)\Big) = \frac{1000}{800-T(x,y,z)}$$
#.
#. Note that there are not thermal expansion effects (i.e. the thermal expansion coefficient is $\alpha=0$).
#. Yet, suprinsingly enough, the problem has analytical solutions for both the temperature and the displacement fields.
#.
#. References:
#. 
#.  * <http://www.code-aster.org/V2/doc/default/fr/man_v/v7/v7.03.100.pdf>
#.  * <https://www.seamplex.com/docs/SP-FI-17-BM-12F2-A.pdf>
#.
#. ## Thermal problem
#.
#. The following input solves the thermal problem over a coarse first-order mesh, writes the resulting
#. temperature distribution into `parallelepiped-temperature.msh`, and prints the $L_2$ error of the numerical
#. result with respect to the analytical solution $T(x,y,z) = 40 - 2x - 3y - 4z$.
#.
#$ $ gmsh -3 parallelepiped.geo -order 1 -clscale 2 -o parallelepiped-coarse.msh
#$ [...]
#$ Info    : 117 nodes 567 elements
#$ Info    : Writing 'parallelepiped-coarse.msh'...
#$ Info    : Done writing 'parallelepiped-coarse.msh'
#$ Info    : Stopped on Fri Dec 10 10:32:30 2021 (From start: Wall 0.0386516s, CPU 0.183052s)
#$ $ feenox parallelepiped-thermal.fee 
#$ 6.18981e-12
#$ $
PROBLEM thermal 3D
READ_MESH parallelepiped-coarse.msh

k = 1     # unitary non-dimensional thermal conductivity

# boundary conditions
BC left    q=+2
BC right   q=-2
BC front   q=+3
BC back    q=-3
BC bottom  q=+4
BC top     q=-4
BC A       T=0

SOLVE_PROBLEM
WRITE_MESH parallelepiped-temperature.msh T

# compute the L-2 norm of the error in the displacement field
Te(x,y,z) = 40 - 2*x - 3*y - 4*z   # analytical solution, "e" means exact
INTEGRATE (T(x,y,z)-Te(x,y,z))^2 RESULT num
PHYSICAL_GROUP bulk DIM 3  # this is just to compute the volume
PRINT num/bulk_volume
#> ![Temperature distribution over the coarse mesh in Gmsh (yes, it is a rainbow pallete)](parallelepiped-temperature.png)
