DifDop Home

DifDop


The Heat Equation

Let us consider the simplest example, the heat equation:

\[ \frac{\partial T}{\partial t}-\frac{\partial^2 T}{\partial x^2}=0; \]

It is found in the directory "examples/heateqn". In this directory, you can find three files that define the problem.

  • The init.xml file where the problem and integration parameters are defined.
    <?xml version="1.0"?>
    <solver nfields="1" naux="1" dt="1e-8" dtout="0.01" tmax="0.1" name="heateqn">
    
      <gsolver 
        method="rkf45"
        step="adaptive"
        eps_abs="1e-6"
        eps_rel="1e-8">
      </gsolver>
    
      <field
        index="1" 
        name="T" 
        init="T_init.h5"> 
      </field>
    
      <mesh 
        x0="0.0" 
        dx="0.005"
        nx="201">
      </mesh>
    
      <params 
        D="1.0">
      </params>
    
    </solver>
    
  • The "init_fields.m" file, which is an octave file that defines the initial condition.
    dx=0.005;
    xx=0.0:dx:1.0;
    T=exp(-(xx-0.5).^2/0.01);
    save('-hdf5',"T_init.h5","T");
    
  • The "heat.c" file where the equation itself is defined.
    #include "mesh.h"
    #include "solver.h"
    #include <stdio.h>
    #include <mxml.h>
    #include <gsl/gsl_errno.h>
    #include <gsl/gsl_odeiv.h>
    
    int fcal_diffn(double t, const double y[], double dydt[], void *params){
      pdsol_rsolver* sv = (pdsol_rsolver *)params;
      pdsol_mesh *rm = sv->rm;
      const double *T;
      double dxp,dxm,D;
      double *DTdt;
      double *S;
      int l;
      T=y;
      DTdt=dydt;
      D=0.1;//sv->params[0].val;
      for(l=1;l<rm->nx-1;l++){
        dxp=rm->xx[l+1]-rm->xx[l];
        dxm=rm->xx[l]-rm->xx[l-1];
        DTdt[l]=((T[l+1]-T[l])/dxp-(T[l]-T[l-1])/dxm)/((dxp+dxm)/2);
      }
      DTdt[0]=0;
      DTdt[rm->nx-1]=0;
      return GSL_SUCCESS;
    }
    
    int main (void){
      double t=0;
      pdsol_rsolver *sv;
      int l;
      sv=pdsol_init_rsolver("init.xml");
      sv->gsv->sys->function=fcal_diffn;
      pdsol_gslrun(sv);
      return 0;
    }
    

If you want to simply run this example by changing some parameters. You can do so by typing:

 octave init_fields.m
 ./heat

if the source code (heat.c) is modified one needs to re-make it by calling:

 make clean;
 make

Fusion Links

Here is a list of links that you might find interesting if you want to learn more about fusion plasmas


This page is currently under construction.