|
User's Heat Equation
Heat Equation
_________
-
α ∇2 U = Ut
- Have a Heat equation to solve? Or any
other math equations? For the first few months of 2009, we are willing
to help you solve them using Calculus-level programming. (Fortran
Calculus is scheduled to be online in mid-2009 so take advantage of
this free offer will it last.) To start, copy and modify the source
code below in a file we'll call it {abc}{123}.fc where {abc} = your
initials and {123} = any number or id; 8 characters max. for filename.
Edit your {abc}{123}.fc file, especially lines starting with a "!"
character.
E-mail
us your {abc}{123}.fc file. We will compile &
execute it and send you the output file. (For more [~60] example
problems to choose from, download FC-Win
program.)
- All "!" characters in columns 1 or 2 must
be deleted before compiling. These "!" were added in order to point out
areas needing work.
- Next, you may modify a copy of this web
page and send it to us for viewing. If accepted, we will post your
webpage showing your problem with solution. If you want people to be
able to contact you, please include your e-mail address on this web
page.
- Please mention our
www.digitalCalculus.com website to others. Thanks!
|
|
|
User's Heat Equation Source Code:
-
For 1-Dimensional Heat Equation use following:
global all
problem HeatPDE
C ------------------------------------------------------------------------
C --- Calculus Programming example: Heat Equation; a PDE Initial
C --- Value Problem solved using Method of Lines.
C ------------------------------------------------------------------------
dynamic U, Ux, Uxx
C
C User parameters ...
! alpha = 2.e-5 ! thermal diffusion coefficienct
! aK = 10. ! W/mK ! thermal conductivity (in the solid)
! h = 25. ! W/m2K ! heat convection at slab surface
! ipoints = 11 ! grid pts. over x-axis
! mpoints = 100 ! grid pts. over t-axis
! tFinal = 6000 ! final time
C
C x-parameter initial settings: x ==> i
! xFinal = 1: xPrint = xFinal/ipoints
C
C t-parameter initial settings: t ==> m
pi = 4*atan(1): mp = mpoints: dt = tFinal/mpoints
allot U( mp), Ux( mp), Uxx( mp)
C
! U0=0
print *, ak, h, U0, dt
C
call xAxis !
end ! Stmt.s not necessary in IVP, but used in BVP versions
model xAxis !
C ... Integrate over x-axis
C
x= 0: xPrt = xPrint: dx = xPrt / 10
Initiate janus; for PDE;
~ equations Uxx/Ux, Ux/U; of x; step dx; to xPrt
do while (x .lt. xFinal)
Integrate PDE; by janus
if( x .ge. xPrt) then
print 79, x, (U( mm), mm = 1, mp)
xPrt = xPrt + xPrint
end if
end do
79 format( 1h , f8.4, /2x, 10(g14.5, 2x))
end
model PDE ! Partial Differential Equation
C ! Method of Lines
! U(1)=(2*h*U0*dx - ak*U(3) + 4*ak*U(2)) / (3*ak + 2*h*dx) ! BC for a slab @ surface
do 20 mm = 2, mpoints ! System of ODEs
tmp = U(mm) - U(mm-1)
Uxx(mm)= 1/alpha * tmp / dt
20 continue
end
-
For 2-Dimensional Heat Equation use following:
global all
problem HeatPDE
C ------------------------------------------------------------------------
C --- Calculus Programming example: Heat Equation; a PDE (2D) Initial
C --- Value Problem solved using Method of Lines.
C ------------------------------------------------------------------------
dynamic U, Ux, Uxx
C
C User parameters ...
! alpha = 2.e-5 ! thermal diffusion coefficienct
! aK = 10. ! W/mK ! thermal conductivity (in the solid)
! h = 25. ! W/m2K ! heat convection at slab surface
! ipoints = 11 ! grid pts. over x-axis
! mpoints = 100 ! grid pts. over t-axis
! tFinal = 6000 ! final time
C
C x-parameter initial settings: x ==> i
! xFinal = 1: xPrint = xFinal/ipoints
C
C y-parameter initial settings: y ==> j
! yFinal = 1: dy = yFinal/(jpoints-1): jp = jpoints
C
C t-parameter initial settings: t ==> m
pi= 4*atan(1): mp = mpoints: dt = tFinal/mpoints
allot U(jp,mp), Ux(jp,mp), Uxx(jp,mp)
C
! U0=0
print *, ak, h, U0, dt
C
call xAxis !
end ! Stmt.s not necessary in IVP, but used in BVP version
model xAxis !
C ... Integrate over x-axis
C settings at t = 0
do 10 mm = 1, mpoints ! Initial Conditions, if any
! xPrt = mm * xPrint: Ux(mm) = 0: Uxx(mm) = 0
! U(mm) = 100.
10 continue
C
x= 0: xPrt = xPrint: dx = xPrt / 10
Initiate janis; for PDE;
~ equations Uxx/Ux, Ux/U; of x; step dx; to xPrt
do while (x .lt. xFinal)
Integrate PDE; by janis
if( x .ge. xPrt) then
do 30 jj = 2, jpoints
y = (jj - 1) * dy
print 79, x, y, (U(jj,mm), mm = 1, mp)
30 continue
xPrt = xPrt + xPrint
end if
end do
79 format( 1h , f8.4, 1x, f8.4, / 2x, 10(g14.5, 1x))
end
model PDE ! Partial Differential Equation
C ! Method of Lines
do 40 jj = 2, jpoints - 1 ! System of ODEs
y = (jj - 1) * dy
do 20 mm = 2, mpoints-1
Ut = (U(jj,mm) - U(jj,mm-1))/dt
Uyy = (U(jj+1,mm) - 2*U(jj,mm) + U(jj-1,mm))/(dy*dy)
Uxx(jj,mm)= 1/alpha * Ut - Uyy
20 continue
! Ux(jj,mp)= ???: Uxx(jj,mp)= ??? ! Initial Conditions, if any
40 continue
end
-
For 3-Dimensional Heat Equation use following:
global all
problem HeatPDE
C ------------------------------------------------------------------------
C --- Calculus Programming example: Heat Equation; a PDE (3D) Initial
C --- Value Problem solved using Method of Lines.
C ------------------------------------------------------------------------
dynamic U, Ux, Uxx
C
C User parameters ...
! alpha = 2.e-5 ! thermal diffusion coefficienct
! aK = 10. ! W/mK ! thermal conductivity (in the solid)
! h = 25. ! W/m2K ! heat convection at slab surface
! ipoints = 11 ! grid pts. over x-axis
! mpoints = 100 ! grid pts. over t-axis
! tFinal = 6000 ! final time
C
C x-parameter initial settings: x ==> i
! xFinal = 1: xPrint = xFinal/ipoints
C
C y-parameter initial settings: y ==> j
! yFinal = 1: dy = yFinal/(jpoints-1): jp = jpoints
C
C z-parameter initial settings: z ==> k
! zFinal = 1: dz = zFinal/(kpoints-1): kp = kpoints
C
C t-parameter initial settings: t ==> m
pi= 4*atan(1): mp = mpoints: dt = tFinal/mpoints
allot U( mp), Ux( mp), Uxx( mp)
C
! U0=0
print *, ak, h, U0, dt
C
call xAxis !
end ! Stmt.s not necessary in IVP, but used in BVP version
model xAxis !
C ... Integrate over x-axis
C
C settings at t = 0
do 10 mm = 1, mpoints ! Initial Conditions, if any
! xPrt = mm * xPrint: Ux(mm) = 0: Uxx(mm) = 0
! U(mm) = 100.
10 continue
C
x= 0: xPrt = xPrint: dx = xPrt / 10
Initiate gemini; for PDE;
~ equations Uxx/Ux, Ux/U; of x; step dx; to xPrt
do while (x .lt. xFinal)
Integrate PDE; by gemini
if( x .ge. xPrt) then
do 30 kk = 2, kpoints
z = (kk - 1) * dz
do 25 jj = 2, jpoints
y = (jj - 1) * dy
print 79, x, y, z, (U(jj,kk,mm), mm = 1, mp)
25 continue
30 continue
xPrt = xPrt + xPrint
end if
end do
79 format( 1h , 3(f8.4, 1x), /2x, 10(g14.5, 1x))
end
model PDE ! Partial Differential Equation
C ! Method of Lines
do 60 kk = 2, kpoints-1 ! System of ODEs
z = (kk - 1) * dz
do 40 jj = 2, jpoints-1
y = (jj - 1) * dy
do 20 mm = 2, mpoints-1
Ut = (U(jj,kk,mm)-U(jj,kk,mm-1))/dt ! approx. partial of U w.r.t. t
Utt = (U(jj,kk,mm+1) - 2*U(jj,kk,mm) + U(jj,kk,mm-1))/dt**2
Uyy = (U(jj+1,kk,mm) - 2*U(jj,kk,mm) + U(jj-1,kk,mm))/dy**2
Uzz = (U(jj,kk+1,mm) - 2*U(jj,kk,mm) + U(jj,kk-1,mm))/dz**2
Uxx(jj,kk,mm)= 1/alpha * Ut - Uyy - Uzz
20 continue
! Ux(jj,kk,mp)= ???: Uxx(jj,kk,mp)= ??? ! Initial Conditions, if any
40 continue
60 continue
end
User's Heat Equation Output:
-
selected output goes here ...
- Visit ODEcalc
for an Ordinary Differential Equations Calculator. Try it and
see the power of Calculus
programming.
HTML code for linking to this page:
- <a
href="http://www.digitalcalculus.com/math-problems/heat-equation.html"><img
align="middle" width="100"
src="http://www.digitalcalculus.com/image/fc-win-icon.gif"/>
<strong>Heat Equation</strong>
</a>; Simulation to Optimization, Tweak Parameters for Optimal
Solution.
Go to top
|
|