🔬
Solver Purpose & Physical Scope
Calculates directional thermal fluxes and rotated effective thermal conductivities for anisotropic materials (e.g., carbon-fiber reinforced composites, wood, and layered PCBs) using 2D tensor rotation formulations: $k(\theta) = k_x \cos^2(\theta) + k_y \sin^2(\theta)$.
📂 Discipline: Conduction
⚡ Precision: IEEE-754 64-bit Real(`real(8)`)
📥 Total Downloads: 247 times
📄 Source File:
anisotropic_cond.f90
📁 calcul/Conduction /
anisotropic_cond.f90
program anisotropic_cond
implicit none
integer :: i
double precision :: kx,ky,kz,Lx,Ly,Lz,Th,Tc,Ac,theta
double precision :: dT,qx,qy,qz,Qx_rate,Qy_rate,Qz_rate,Rx,Ry,Rz
double precision :: keff,th_rad,ths,keffs,pi
pi=3.14159265358979d0
read(*,*) kx; read(*,*) ky; read(*,*) kz
read(*,*) Lx; read(*,*) Ly; read(*,*) Lz
read(*,*) Th; read(*,*) Tc; read(*,*) Ac; read(*,*) theta
dT=abs(Th-Tc)
th_rad=theta*pi/180d0
qx=kx*dT/Lx; qy=ky*dT/Ly; qz=kz*dT/Lz
Qx_rate=qx*Ac; Qy_rate=qy*Ac; Qz_rate=qz*Ac
Rx=Lx/(kx*Ac); Ry=Ly/(ky*Ac); Rz=Lz/(kz*Ac)
keff=kx*cos(th_rad)**2+ky*sin(th_rad)**2
write(*,'(A)') '============================================'
write(*,'(A)') ' ANISOTROPIC THERMAL CONDUCTION'
write(*,'(A)') '============================================'
write(*,'(A)') ''
write(*,'(A)') '--- INPUTS ---'
write(*,'(A,F12.4,A)') ' k_x = ',kx,' W/mK'
write(*,'(A,F12.4,A)') ' k_y = ',ky,' W/mK'
write(*,'(A,F12.4,A)') ' k_z = ',kz,' W/mK'
write(*,'(A,F10.4,A)') ' L_x = ',Lx,' m'
write(*,'(A,F10.4,A)') ' L_y = ',Ly,' m'
write(*,'(A,F10.4,A)') ' L_z = ',Lz,' m'
write(*,'(A,F10.2,A)') ' T_hot = ',Th,' C'
write(*,'(A,F10.2,A)') ' T_cold = ',Tc,' C'
write(*,'(A,F10.4,A)') ' Cross Section A = ',Ac,' m2'
write(*,'(A,F10.2,A)') ' Rotation theta = ',theta,' deg'
write(*,'(A)') ''
write(*,'(A)') '--- DIRECTIONAL PROPERTIES ---'
write(*,'(A,F12.4,A)') ' Heat flux q_x = ',qx,' W/m2'
write(*,'(A,F12.4,A)') ' Heat flux q_y = ',qy,' W/m2'
write(*,'(A,F12.4,A)') ' Heat flux q_z = ',qz,' W/m2'
write(*,'(A,F12.4,A)') ' Heat rate Q_x = ',Qx_rate,' W'
write(*,'(A,F12.4,A)') ' Heat rate Q_y = ',Qy_rate,' W'
write(*,'(A,F12.4,A)') ' Heat rate Q_z = ',Qz_rate,' W'
write(*,'(A,F12.6,A)') ' Resistance R_x = ',Rx,' K/W'
write(*,'(A,F12.6,A)') ' Resistance R_y = ',Ry,' K/W'
write(*,'(A,F12.6,A)') ' Resistance R_z = ',Rz,' K/W'
write(*,'(A)') ''
write(*,'(A)') '--- ANISOTROPY RATIOS ---'
write(*,'(A,F10.4)') ' k_x / k_y = ',kx/ky
write(*,'(A,F10.4)') ' k_x / k_z = ',kx/kz
write(*,'(A,F10.4)') ' k_y / k_z = ',ky/kz
write(*,'(A,F12.4,A)') ' k_eff(theta) in x-y = ',keff,' W/mK'
write(*,'(A)') ''
write(*,'(A)') '--- ROTATED k_eff (x-y PLANE) ---'
write(*,'(A)') ' Tensor rotation: k_eff = kx*cos^2(theta) + ky*sin^2(theta)'
write(*,'(A)') ''
write(*,'(A)') '--- ANGLE SWEEP ---'
write(*,'(A)') ' theta[deg] k_eff[W/mK] q[W/m2] Q[W]'
write(*,'(A)') ' -----------------------------------------------'
do i=1,19
ths=dble(i-1)*5d0
keffs=kx*cos(ths*pi/180d0)**2+ky*sin(ths*pi/180d0)**2
write(*,'(2X,F8.1,4X,F10.4,4X,F10.4,4X,F10.4)') ths,keffs,keffs*dT/Lx,keffs*dT/Lx*Ac
enddo
write(*,'(A)') ''
write(*,'(A)') '--- CORRELATIONS ---'
write(*,'(A)') ' k_eff(theta) = kx*cos^2(theta) + ky*sin^2(theta)'
write(*,'(A)') ' Off-diagonal: k_xy = (kx-ky)*sin(theta)*cos(theta)'
write(*,'(A)') ' 3D rotation requires full tensor transformation'
write(*,'(A)') ' Ref: Incropera Appendix A, Ozisik Ch.1'
write(*,'(A)') ' Chung, Thermal Conductivity of Composites'
end program anisotropic_cond
💻 How to Compile & Run Locally
1. Compilation (GNU Fortran / Intel oneAPI):
gfortran -O3 anisotropic_cond.f90 -o anisotropic_cond
2. Execution with input.txt redirection:
anisotropic_cond < input.txt
📄 Sample input.txt File Structure
Sample Data:
7.0 0.8 0.8 0.01 0.01 0.01 100.0 25.0 0.001 0.0
Parameter Description:
Conductivity k_x [W/mK] Conductivity k_y [W/mK] Conductivity k_z [W/mK] Thickness L_x [m] Thickness L_y [m] Thickness L_z [m] Hot surface temperature T_hot [C] Cold surface temperature T_cold [C] Cross-section area A [m2] Rotation angle theta [deg]