Commit fef656df authored by Christian Engwer's avatar Christian Engwer

Added slightly more sophisticated vtkout function that is also used for

 parallel runs of the finite volume example. The writePVD shell script makes a
 paraview video data file which is documented in the howto. The script is copied
 from the dune-femhowto.

creadits to Martin Drohmann

[[Imported from SVN: r235]]
parent 342705e5
......@@ -48,12 +48,12 @@ void timeloop (G& grid, double tend, int lmin, int lmax)
for (int i=grid.maxLevel(); i<lmax; i++)
{
if (grid.maxLevel()>=lmax) break;
finitevolumeadapt(grid,mapper,c,lmin,lmax,0);
finitevolumeadapt(grid,mapper,c,lmin,lmax,0); /*@\label{afv:in}@*/
initialize(grid,mapper,c);
}
// write initial data
vtkout(grid,c,"concentration",0);
vtkout(grid,c,"concentration",0,0);
// variables for time, timestep etc.
double dt, t=0;
......@@ -78,7 +78,7 @@ void timeloop (G& grid, double tend, int lmin, int lmax)
if (t >= saveStep)
{
// write data
vtkout(grid,c,"concentration",counter);
vtkout(grid,c,"concentration",counter,t);
// increase counter and saveStep for next interval
saveStep += saveInterval;
......@@ -86,17 +86,20 @@ void timeloop (G& grid, double tend, int lmin, int lmax)
}
// print info about time, timestep size and counter
std::cout << "s=" << grid.size(0) << " k=" << k << " t=" << t << " dt=" << dt << std::endl;
std::cout << "s=" << grid.size(0)
<< " k=" << k << " t=" << t << " dt=" << dt << std::endl;
// for unstructured grids call adaptation algorithm
if( Dune :: Capabilities :: IsUnstructured<G> :: v )
{
finitevolumeadapt(grid,mapper,c,lmin,lmax,k);
finitevolumeadapt(grid,mapper,c,lmin,lmax,k); /*@\label{afv:ad}@*/
}
}
// write last time step
vtkout(grid,c,"concentration",counter);
vtkout(grid,c,"concentration",counter,tend);
// write
}
//===============================================================
......
......@@ -43,14 +43,14 @@ void timeloop (const G& grid, double tend)
// initialize concentration with initial values
initialize(grid,mapper,c); /*@\label{fvc:init}@*/
vtkout(grid,c,"concentration",0);
vtkout(grid,c,"concentration",0,0.0);
// now do the time steps
double t=0,dt;
int k=0;
const double saveInterval = 0.1;
double saveStep = 0.1;
int counter = 0;
int counter = 1;
while (t<tend) /*@\label{fvc:loop0}@*/
{
......@@ -67,7 +67,7 @@ void timeloop (const G& grid, double tend)
if (t >= saveStep)
{
// write data
vtkout(grid,c,"concentration",counter);
vtkout(grid,c,"concentration",counter,t);
// increase counter and saveStep for next interval
saveStep += saveInterval;
......@@ -75,11 +75,12 @@ void timeloop (const G& grid, double tend)
}
// print info about time, timestep size and counter
std::cout << "s=" << grid.size(0) << " k=" << k << " t=" << t << " dt=" << dt << std::endl;
std::cout << "s=" << grid.size(0)
<< " k=" << k << " t=" << t << " dt=" << dt << std::endl;
} /*@\label{fvc:loop1}@*/
// output results
vtkout(grid,c,"concentration",counter); /*@\label{fvc:file}@*/
vtkout(grid,c,"concentration",counter,tend); /*@\label{fvc:file}@*/
}
//===============================================================
......@@ -95,7 +96,7 @@ int main (int argc , char ** argv)
try {
using namespace Dune;
// use unitcube from grids
// use unitcube from dgf grids
std::stringstream dgfFileName;
dgfFileName << "grids/unitcube" << GridType :: dimension << ".dgf";
......
......@@ -42,21 +42,41 @@ void partimeloop (const G& grid, double tend)
// initialize concentration with initial values
initialize(grid,mapper,c);
vtkout(grid,c,"pconc",0);
vtkout(grid,c,"pconc",0,0.0,grid.comm().rank());
// now do the time steps
double t=0,dt;
int k=0;
const double saveInterval = 0.1;
double saveStep = 0.1;
int counter = 1;
while (t<tend)
{
// augment time step counter
k++;
// apply finite volume scheme
parevolve(grid,mapper,c,t,dt);
// augment time
t += dt;
// check if data should be written
if (t >= saveStep)
{
// write data
vtkout(grid,c,"pconc",counter,t,grid.comm().rank());
//increase counter and saveStep for next interval
saveStep += saveInterval;
++counter;
}
// print info about time, timestep size and counter
if (grid.comm().rank()==0) /*@\label{pfc:rank0}@*/
std::cout << "k=" << k << " t=" << t << " dt=" << dt << std::endl;
if (k%20==0) vtkout(grid,c,"pconc",k/20);
}
vtkout(grid,c,"pconc",k/20);
vtkout(grid,c,"pconc",counter,tend,grid.comm().rank());
}
//===============================================================
......
......@@ -4,11 +4,20 @@
#include <stdio.h>
template<class G, class V>
void vtkout (const G& grid, const V& c, const char* name, int k)
void vtkout (const G& grid, const V& c, const char* name, int k, double time=0.0, int rank=0)
{
Dune::VTKWriter<typename G::LeafGridView> vtkwriter(grid.leafView());
char fname[128];
char sername[128];
sprintf(fname,"%s-%05d",name,k);
sprintf(sername,"%s.series",name);
vtkwriter.addCellData(c,"celldata");
vtkwriter.write(fname,Dune::VTKOptions::ascii);
if ( rank == 0)
{
std::ofstream serstream(sername, (k==0 ? std::ios_base::out : std::ios_base::app));
serstream << k << " " << fname << ".vtu " << time << std::endl;
serstream.close();
}
}
#! /bin/bash
# read parameters
if [ -z "$1" ]; then
echo "Usage: writePVD [file_prefix]"
exit
fi
param=$1
SERIESFILE="${param##*/}"
DIR="${param%$SERIESFILE}"
DIR="${DIR:-.}"
FILEPREFIX="${SERIESFILE%.*}"
SERIESFILE=$1
echo $SERIESFILE $DIR $FILEPREFIX
THISDIR=$PWD
# cd $DIR
# output file name
OUTPUTFILE=$FILEPREFIX.pvd
#write header
echo "Writing file: $OUTPUTFILE"
echo "<?xml version=\"1.0\"?>" > $OUTPUTFILE
echo "<VTKFile type=\"Collection\" version=\"0.1\" byte_order=\"LittleEndian\">" >> $OUTPUTFILE
echo " <Collection>" >> $OUTPUTFILE
echo "" >> $OUTPUTFILE
# Set loop separator to end of line
BAKIFS=$IFS
IFS=$(echo -en "\n\b")
# the following applied to the series file
exec 3<&0
exec 0<$SERIESFILE
while read line
do
# file name and time stamp
NAME=$(echo $line | cut -d " " -f 2)
TIME=$(echo $line | cut -d " " -f 3)
# strip unnecessary path from file name
while ! test -e $NAME
do
NAME="${NAME#*/}"
if [ $NAME == "" ];
then # file does not exists
NAME=$(echo $line | cut -d " " -f 2)
echo $NAME " not found!"
exit
fi
done
# write line into output file
echo "<DataSet timestep=\"$TIME\" group=\"\" part=\"0\" file=\"$NAME\"/>" >> $OUTPUTFILE
done
exec 0<&3
# restore $IFS which was used to determine what the field separators are
IFS=$BAKIFS
# write the end of file
echo "" >> $OUTPUTFILE
echo " </Collection>" >> $OUTPUTFILE
echo "</VTKFile>" >> $OUTPUTFILE
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment