This example demonstrates how to create current streamlines
with in a torso mesh (BE) using domain solutions. This makes use of
a Perl routine in the command file. In the main example directory there is an example_cmiss
command file which can be run using 'cmiss' and will do both the back and front
end calculations/visualisation. There is also an example_cmgui command file which
performs the visualisation in cmgui only.
if( !$example ) {
$example = ".";
}
$FBnone = "\033[0m";
$Fgreen = "\033[32m";
$Fblue = "\033[34m";
print( "$Fgreen" );
# This example is designed to go through the whole process of getting current lines
# within a boundary element torso using streamlines. This is achieved by putting the torso
# inside a box and then defining a grid over the box. A domain solution is then performed
# at each grid point to get an element based representation of the potential gradient fields.
# A Perl routine is then used to combine the x, y and z components of the gradient field into
# a single exelem file.
print( "$FBnone" );
print( "$Fgreen" );
# Set up the three regions for the example
print( "$FBnone" );
$Torso = 1;
$Heart = 2;
$Box = 3;
print( "$Fgreen" );
# Set up the two classes for the example
print( "$FBnone" );
$Dipole = 1;
$Grid = 2;
print( "$Fgreen" );
# Define the problem setup - 3d with 3 regions
print( "$FBnone" );
fem define para;r;$example/torsotimedep
fem define regi;r;$example/torsouno
fem define coor 3,1
print( "$Fgreen" );
# Read in the geometry of the 3 regions
print( "$FBnone" );
fem define node;r;$example/torsotimedep region $Torso,$Heart
fem define node;r;$example/volume reg $Box
fem define base;r;$example/volume
fem define elem;r;$example/torsotimedep region $Torso,$Heart
fem define elem;r;$example/volume reg $Box
print( "$Fgreen" );
# Define a finite difference grid over the box that encloses the
# torso. The potential field will be evaluated at these points.
print( "$FBnone" );
fem define grid;r;$example/45k reg $Box
fem update grid geom reg $Box
print( "$Fgreen" );
# Define a quasistatic Laplace equation in the torso regions
# along with the material parameters and initial conditions etc.
print( "$FBnone" );
fem define equa;r;$example/torsotimedep region $Torso,$Heart class $Dipole lock
fem define mate;r;$example/torsotimedep region $Torso,$Heart
fem define init;r;$example/torsotimedep region $Torso,$Heart
fem define sour;r;$example/test region $Heart
fem define coup;r;$example/torsotimedep
fem define solv;r;$example/torsotimedep coupled region $Torso,$Heart
print( "$Fgreen" );
# A grid equation needs to be defined so we have a grid class in
# which to store the potentials. Here a Laplace equation is used for
# this purpose.
print( "$FBnone" );
fem define equa;r;$example/volume cl $Grid reg $Box
print( "$Fgreen" );
# Solve the torso problem. In this case it has been set up to solve for
# only the time step that we want to export the potentials from.
print( "$FBnone" );
fem solve coupled region $Torso,$Heart class $Dipole
print( "$Fgreen" );
# Export the geometry of the torso and epicardial surfaces
print( "$FBnone" );
fem export node;hearttorso as hearttorso
fem export elem;hearttorso as hearttorso
print( "$Fgreen" );
# Export the potential field over the torso and epicardial surfaces
print( "$FBnone" );
fem export node;field as hearttorso field
fem export elem;field as hearttorso field
print( "$Fgreen" );
# Perform a boundary element domain solution at each of the grid points
# to find the value of the potential gradient fields. There are two easy ways to get
# better, smoother streamlines but both take substantially longer.
# They involve:
# Increasing grid density
# Using adaptive integration
# This evaluation provides linear speedup with accurate answers when
# mulitprocessing (if not it certainly has done in the past).
print( "$FBnone" );
fem evaluate solution grid element 345 optimise iy 1 gradient
print( "$Fgreen" );
# Export the geometry of the box that encloses the torso
print( "$FBnone" );
fem export node;geom as test reg $Box
fem export elem;geom as test reg $Box
print( "$Fgreen" );
# Export the element based fields that contains the potential gradient fields evaluated
# at all grid points in the torso. All points outside the torso have been set to
# 0.0mV based of the C(P) coefficient that is calculated by the boundary element
# method.
print( "$FBnone" );
fem export elem;xpote as test field reg $Box class $Grid iy 2
fem export elem;ypote as test field reg $Box class $Grid iy 3
fem export elem;zpote as test field reg $Box class $Grid iy 4
print( "$Fblue" );
sub WriteGradPhiExelem {
# This Perl subroutine is used to combine three potentential gradient
# exelem files into a single three component field within a single
# exelem file. This routine requires filenames for the three input
# exelems to be passed and also a filename for the output exelem.
( $DXFILENAME, $DYFILENAME, $DZFILENAME, $OUTPUTNAME ) = @_;
my $Found = 0;
my $Found2 = 0;
my $First = 1;
my $InputLine = ' ';
CORE::open( OUTPUTFILE,">$OUTPUTNAME" ) or die( "ERROR: Unable to open output file\n" );
CORE::open( INPUTFILE,"<$DXFILENAME" ) or die( "ERROR: Unable to open DX input file\n" );
while( <INPUTFILE> ) {
if( /Nodes:/ ) {
$Inputline = $_;
$Found2 = 0;
CORE::open( INPUTFILE2,"<$DYFILENAME" ) or die( "ERROR: Unable to open DY input file\n" );
while( <INPUTFILE2> ) {
if( /Nodes:/ ) {
$Found2 = 0;
}
if( $Found2 ) {
print( OUTPUTFILE "$_" );
}
if( /Values:/ ) {
$Found2 = 1;
}
}
CORE::close( INPUTFILE2 );
$Found2 = 0;
CORE::open( INPUTFILE2,"<$DZFILENAME" ) or die( "ERROR: Unable to open DZ input file\n" );
while( <INPUTFILE2> ) {
if( /Nodes:/ ) {
$Found2 = 0;
}
if( $Found2 ) {
print( OUTPUTFILE "$_" );
}
if( /Values:/ ) {
$Found2 = 1;
}
}
CORE::close( INPUTFILE2 );
print( OUTPUTFILE "$Inputline" );
}
if( !$Found ) {
print( OUTPUTFILE $_ );
} elsif( $First ) {
print( OUTPUTFILE ' 1) gradphi, field, rectangular cartesian, #Components=3'."\n");
print( OUTPUTFILE ' x. l.Lagrange*l.Lagrange*l.Lagrange, no modify, grid based.'."\n");
$First = 0;
}
if( /Fields/ ) {
$Found = 1;
}
if( /xi1/ ) {
print( OUTPUTFILE $_ );
print( OUTPUTFILE ' y. l.Lagrange*l.Lagrange*l.Lagrange, no modify, grid based.'."\n");
print( OUTPUTFILE $_ );
print( OUTPUTFILE ' z. l.Lagrange*l.Lagrange*l.Lagrange, no modify, grid based.'."\n");
print( OUTPUTFILE $_ );
}
if(( /Element/ ) && $Found ) {
print( OUTPUTFILE $_ );
$Found = 0;
}
}
CORE::close( INPUTFILE );
CORE::close( OUTPUTFILE );
} # WriteGradPhiExelem
print( "$FBnone" );
print( "$Fgreen" );
# Set up the filenames to pass to the Perl subroutine
print( "$FBnone" );
$DXFILENAME = "xpote.exelem";
$DYFILENAME = "ypote.exelem";
$DZFILENAME = "zpote.exelem";
$OUTPUTNAME = "GradPhi.exelem";
print( "$Fgreen" );
# Call the Perl subroutine to generate the exelem file that can be read
# into CMGUI to create the streamlines.
print( "$FBnone" );
&WriteGradPhiExelem( $DXFILENAME, $DYFILENAME, $DZFILENAME, $OUTPUTNAME );