{  GOS_FlexPDE_2D_Beam_01.pde
}

title 'Beam 3D heat conduction'

select
    regrid=off { use fixed grid }
    ngrid=30     { smaller grid for quicker run }

coordinates
    cartesian2

variables
    Tp(threshold=0.01)   { the temperature variable, with approximate size  }
    Tpc
	TpcSurface
    TpcSurfaceMax

time 0 to 100000 by 60 { establish time range and initial timestep }

definitions

    timeMax = 50400
    Ta = 20 + 273
    bottomTemp = Ta
    surfaceTemp = Ta
    length = 1					{ length }
    wide = 1.2					{ width }
    thick = 0.2					{ thickness }
    beamTopY = thick/2
    beamTopLowerY = beamTopY - 0.05

    K = 1.3  						{ Thermal conductivity  }
    H = 6
	Cp = 880 
    m = 2400 ! this is really rho, the density of the material !
    C = Cp * m

 Q = 20 * Cp * Sin(2*pi* t / 36000)  * upulse(t , t - 18000)   { Thermal source }

! Cannot get a max along this X axis
TpcSurfaceFunc = Passive( IF ((y <= beamTopY) AND (y >= beamTopLowerY)) THEN Tpc ELSE TpcSurfaceMax )
! Can get a max of Tpc, but this is not what I want
! TpcSurfaceFunc(y1) = Tpc
TpcSurfaceMaxFunc = Passive( max(TpcSurface, TpcSurfaceMax) )


initial values	
    Tp = 20 + 273
    Tpc = 20

	TpcSurface = 20
	TpcSurfaceMax = 20

equations
    Tp: div(k*grad(Tp)) + Q = C * dt(Tp) 	{ the heat equation }  
then
	Tpc: Tpc = Passive( Tp - 273 )
then
	TpcSurface: TpcSurface = TpcSurfaceFunc
then
	TpcSurfaceMax: TpcSurfaceMax = TpcSurfaceMaxFunc
  
 boundaries

    Region 'BeamRegion'    							{ Define full domain boundary }
       start "Beam" (-wide/2, -thick/2)
		 	natural(Tp) = -H * 0.5 *(Tp - Ta)  		{ Set convection value for Bottom }
         line to  (wide/2, -thick/2)  			
		 	natural(Tp) = -H * 1 * (Tp - Ta)  		{ Set convection value for Right }
         line to  (wide/2, thick/2)
		 	natural(Tp) = -H * 1 *(Tp - Ta)  		{ Set convection value for Top }
         line to (-wide/2 , thick/2)
		 	natural(Tp) = -H * 1 * (Tp - Ta)  		{ Set convection value for Left }
         line to (-wide/2, -thick/2)


monitors
  for cycle=1
    contour(Tpc) on "BeamRegion" as "Beam Temp" ! painted ! range=(0,tmax)
report TpcSurface as "TpcSurface"   
 !   elevation(Tpc) from (-wide,0,0) to (wide,0,0) as "X-Axis Temp" ! range=(0,tmax)

plots
  for t = timeMax 
    contour(Tpc) on "BeamRegion" as "Beam Temp at Tmax" painted ! range=(0,tmax)

    elevation(Tpc) from (-wide/2, thick/2) to (wide/2, thick/2) as "Top Surface Temp"


histories
    history(Tpc) at 
                    (0, thick/2)  ! range=(0,tmax)
                    (0, 0)  ! range=(0,tmax)
                    (0, -thick/2)  ! range=(0,tmax)

report TpcSurface as "TpcSurface"   
report TpcSurfaceMax as "TpcSurfaceMax"   



end