TITLE 'PISA simulation' COORDINATES cartesian3 VARIABLES temp voltage SELECT { method controls } ! threads = 2 ! multithreading, sometimes not helpful notify_done = on prefer_stability = on prefer_speed = off upwind = on !errlim = 0.1 DEFINITIONS ! o stands for outer, i for innerd everything in SI units !!GEOMETRY!! Ro=7.6e-3 ! outer radius Ri=6e-3 ! inner radius s = (Ro^2 - Ri^2) * PI lencyl= 70e-3 ! length of cylinder hole_bottom=lencyl / 2 - 0.75e-3 !0.4 ! lower coordiante of the hole hole_top=lencyl / 2 + 0.75e-3 !0.6 ! upper coordinate of the hole, hole_top minus hole_bottom is height of hole opening_angle=4.776 ! opening angle !!PHYSICAL PROPERTIES!! temp0=300 ko=54 ! thermal conductivity of metal ki=1e-32 ! thermal conductivity of air or filling rhoo=8.57e6 ! density of material rhoi=1e-32 sigmao=1 / (6.58e6) * temp0 / abs(temp) !* pi * Ro^2! electrical conductivity of metal ! Sheet resistance??? sigmai=1e-32 cpo=0.265 ! specific heat capacity of metal cpi=1e-32 Q=0 ! Heat source or sink !!FURTHER PARAMETERS!! k=ki ! standard thermal conductivity in case that nothing is specified rho=rhoi cp=cpi sigma = sigmai voltage_max = 20 * ramp(t - 0.2, 1, 0, 1e-32) !* ramp(t-0.01, 0, 1, 1e-32) !* (z/lencyl) voltage_up = voltage_max / 2 voltage_down = - voltage_max / 2 ! voltage = (20 * ramp(t - 0.02, 1, 0, 1e-32)) * (z/lencyl) !Transfer mesh !transfermesh('transferm3D_1.xfr') INITIAL VALUES temp=500 !voltage = voltage_down + voltage_max * (z/lencyl) !* ramp(t - 0.02, 1, 0, 1e-32) EQUATIONS voltage: div( -sigma * grad(voltage)) = 0 !temp: div( -k * grad(temp)) - sigma * (magnitude( -grad(voltage))^2) / s^2 - 0.8 * 5.67e-8 * (2 * Ro - 2 * Ri)* (temp^4 - 300^4) / s = dt(temp) * cp * rho ! time dependence with dt(temp)... temp: div(-k * grad(temp)) - sigma * (magnitude( -grad(voltage))^2) / lencyl - 0.8 * 5.67e-8 * (2 * Ro - 2 * Ri)* (temp^4 - 300^4) / s = dt(temp) * cp * rho ! CONSTRAINTS { Integral constraints } EXTRUSION surface "cyl_bottom" z=0 layer "lower" surface "hole_bottom" z=hole_bottom layer "middle" surface "hole_top" z=hole_top layer "upper" surface "cyl_top" z=lencyl BOUNDARIES SURFACE "cyl_bottom" value(temp)=temp0 value(voltage)=voltage_down SURFACE "cyl_top" value(temp)=temp0 value(voltage)=voltage_up REGION 1 layer "lower" k=ko cp=cpo rho=rhoo sigma=sigmao layer "upper" k=ko cp=cpo rho=rhoo sigma=sigmao layer "middle" k=ko cp=cpo rho=rhoo sigma=sigmao START (0,Ro) natural(temp)=temp0 !natural(voltage)=0 arc(center=0,0) angle=360 to close REGION 2 layer "lower" k=ko cp=cpo rho=rhoo sigma=sigmao layer "upper" k=ko cp=cpo rho=rhoo sigma=sigmao layer "middle" k=ko cp=cpo rho=rhoo sigma=sigmao START (-SIN(opening_angle/2 DEGREES)*Ro,-COS(opening_angle/2 DEGREES)*Ro) arc(center=0,0) angle=opening_angle line to (SIN(opening_angle/2 DEGREES)*Ri,-COS(opening_angle/2 DEGREES)*Ri) arc(center=0,0) angle=-opening_angle natural(temp)=temp0 line to close REGION 3 layer "lower" k=ki cp=cpi rho=rhoi sigma=sigmai layer "middle" k=ki cp=cpi rho=rhoi sigma=sigmai layer "upper" k=ki cp=cpi rho=rhoi sigma=sigmai START (0,Ri) natural(temp)=temp0 !natural(voltage)=voltage_max arc(center=0,0) angle=360 to close time 0 to 1 by 0.1! needed for transient problems, Plots must be modified so that a specific time is given MONITORS for cycle=10 contour(temp) on x=0 contour(voltage) on x=Ri SUMMARY !elevation(voltage) from t report(voltage) PLOTS for cycle=10 for t = 0 by 1to 0.1 ! contour(temp, t) on x=0 as "temp on x=0" painted ! contour(voltage, t) on x=Ri as "voltage on x=Ri" painted contour(temp, t) on z=0.5 as "temp on z=0.5" painted PNG ! contour(k) on x=0 as "k on x=0" painted ! contour(k*dz(temp)) on x=0 painted ! contour(k) on z=0.5 as "k on z=0.5" painted ! elevation(temp) from (0,-Ro,0) to (0,-Ro,lencyl) as "temp on side with hole" ! elevation(temp) from (0,Ro,0) to (0,Ro,lencyl) as "temp on backside" contour(voltage) on x=Ri painted PNG transfer(temp) file="transferm3D.xfr" HISTORIES history(temp) at (Ri, 0, lencyl/2) as "history for temp on x=0" PNG history(voltage) at (Ri, 0, lencyl/2) as "history of voltage on x=Ri" PNG END