ff3d-users
[Top][All Lists]

## [ff3d-users] Could you help me confirm a time-dependent PDE?

 From: Robert Li Subject: [ff3d-users] Could you help me confirm a time-dependent PDE? Date: Tue, 24 Feb 2004 15:20:15 -0800 (PST)

```Dear Stephane Del Pino,

I tested ff3d with a realistic thermal PDE. However,
I found that the result of ff3d was different from
that of SVHEAT and TEMP/W, which are also FEM solution
software. But the other two have close result, or
almost same.

I am not sure about how to write PDE in input file.
For example, It seems that -div((f2)*grad(T)) = 0; and
-div((f2)*grad(T)) = 0; are different. But I don't
know why?

Could you help me confirm the following
time-dependent PDE?

1.Real Heat flow Equation
dx(Kx*dx(T)) + dy(Ky*dy(T)) +
dz(Kz*dz(T))=(Cp+Lf*Theat*Theatu1)*(dt(T))
where, T=Temperature,
dt=time

2.I expressed it in ff3d input file as follows:

where,
Tn = previously-computed temperature
dt = time step 1 day
Lf=Latent heat of fusion of water
K=conductivity
Cp=volumetric heat capacity
Theatu1=slope of unfrozen water content
Theat=volumetric water content
function        factor=((dt)/((Cp)+(Lf)*(Theat)*(Theatu1)));
function        f2=(K*factor);

Is it correct?

Thanks!

Robert

I used some selfdefind functions in the following
inputfile for conductivity, etc.

>>>>>>>>>>>input file>>>>>>>>

femfunction     Tn(M) = initemperature();                               //T0;
when
t=0; T
function                Told = 0;                       //T0; when t=0; T
double                  i = 1;                          //loop times
double                  dt = (1);                       //time step 1 day
double                  Lf=334000.0;
//Latent heat of fusion of
water
function        K=(conductivity(Tn,-1) );                       //conductivity
in z direction
function        Cp=(vheatcapacity(Tn,-1));                      //volumetric
heat capacity
function        Theatu1=(slopeofunfwc(Tn, Told, -1));
//use new slope method
function        Theat=theat();

function        factor=((dt)/((Cp)+(Lf)*(Theat)*(Theatu1)));
function        f2=(K*factor);

//Get initial value T0, when time t=0;
solve(T) in  M cg(epsilon=1E-6){
pde(T)

-dx((K)*dx(T))-dy((K)*dy(T))-dz((K)*dz(T))  = 0;

T = -10 on  2;
T = 2   on  1;

}
//Save result
save (medit , ".\\temp\\Model9xxxx\\model.000" , T, M,
dos ) ;
save (medit , ".\\temp\\Model9xxxx\\model.000",  M,
dos ) ;

Tn= T;

//solve
do {
cout<<"I-----------------"<<i<<"\n";
solve(T) in  M cg(epsilon=1E-6){
pde(T)

T = boundarytemperature(i) * correctvalue(Tn) on 2;
T = 2 on 1;

}

//Save result
save (medit , ".\\temp\\Model9xxxx\\model.00".i , T,
M, dos ) ;
save (medit , ".\\temp\\Model9xxxx\\model.00".i,  M,
dos ) ;

i=i+1;

// correct value here
Told = Tn;
Tn = T;

} while (i<=365);

__________________________________
Do you Yahoo!?
Yahoo! Mail SpamGuard - Read only the mail you want.
http://antispam.yahoo.com/tools```
```//load initial parameters

femfunction     Tn(M) = initemperature();                               //T0;
when t=0; T
function                Told = 0;                       //T0; when t=0; T
double                  i = 1;                          //loop times
double                  dt = (1);                       //time step 1 day
double                  Lf=334000.0;
//Latent heat of fusion of water
function        K=(conductivity(Tn,-1) );                       //conductivity
in z direction
function        Cp=(vheatcapacity(Tn,-1));                      //volumetric
heat capacity
function        Theatu1=(slopeofunfwc(Tn, Told, -1));   //use new slope method
function        Theat=theat();

function        factor=((dt)/((Cp)+(Lf)*(Theat)*(Theatu1)));
function        f2=(K*factor);

//Get initial value T0, when time t=0;
solve(T) in  M cg(epsilon=1E-6){
pde(T)

-dx((K)*dx(T))-dy((K)*dy(T))-dz((K)*dz(T))  = 0;

T = -10 on  2;
T = 2   on  1;

}
//Save result
save (medit , ".\\temp\\Model9xxxx\\model.000" , T, M, dos ) ;
save (medit , ".\\temp\\Model9xxxx\\model.000",  M, dos ) ;

Tn= T;

//solve
do {
cout<<"I-----------------"<<i<<"\n";
solve(T) in  M cg(epsilon=1E-6){
pde(T)

T = boundarytemperature(i) * correctvalue(Tn) on 2;
T = 2 on 1;

}

//Save result
save (medit , ".\\temp\\Model9xxxx\\model.00".i , T, M, dos ) ;
save (medit , ".\\temp\\Model9xxxx\\model.00".i,  M, dos ) ;

i=i+1;

// correct value here
Told = Tn;
Tn = T;

} while (i<=365);```