Output files produced by Montjoie

In order to visualize outputs produced by Montjoie, it is possible to use Matlab, Python, Medit or Paraview. Below, we have listed for each type of output, which software can be used :

SismoGrille Python, Matlab, Paraview
SismoGrille3D Python, Matlab, Paraview
SismoPlane Python, Matlab, Paraview
SismoLine Python, Matlab, Paraview
SismoCircle Python, Matlab, Paraview
SismoCircleAxi Python, Matlab, Paraview
SismoPointsFile Python, Matlab, Paraview
SismoPoint Python, Matlab
SismoMeshSurfacic Medit, Paraview
SismoMeshVolumic Medit, Paraview
FileRcs Python, Matlab

Matlab/Python

Matlab/Python scripts are placed in directory Matlab (Python functions are actually regrouped in the single file visuND.py) :

For example, if you want to plot a solution with Matlab, you will type :

cd MATLAB
[X, Y, Z, coor, V] = loadND('../totalCarre_U0.dat');
plot2dinst(X, Y, real(V));

By using ipython -pylab, the commands are similar :

ipython -pylab
from visuND import *
[X, Y, Z, coor, V] = loadND('../totalCarre_U0.dat');
plot2dinst(X, Y, real(V));

For outputs produced by specifying SismoGrille, SismoGrille3D, SismoPlane, SismoLine, SismoCircle, SismoCircleAxi, SismoPointsFile, the function loadND has to be considered, whereas in the case of SismoPoint (for time-domain simulations) or FileRcs, the file can be read directly with the command load (loadtxt in Python).

Paraview

If you wish to visualize results with Paraview, don't forget to insert in the data file one of these lines

OutputFormat = VTK FLOAT
OutputFormat = ASCII_VTK FLOAT
OutputFormat = VTK DOUBLE
OutputFormat = ASCII_VTK DOUBLE

Montjoie should produce Vtk files (extension .vtk), that you can read in paraview (by clicking on File, Open). For 1-D outputs (on lines, circles), you can apply the filter "Plot Data" to see the 1-D curve. If you specify VTK after OutputFormat, the datas will be written in binary, and will take less memory. By specifying float, the datas are written in single precision (less memory is used).

Classes for outputs in Montjoie

Public methods of GridInterpolation

GetElementNumber returns the element where the point i is
GetLocalCoordinate returns the local coordinates of point i in the reference element
GetGlobalCoordinate returns the global coordinates of point i
GetTheta returns the parameter theta for the i-th point (for axisymmetric outputs)
GetDFjm1 returns the inverse of the jacobian matrix for the i-th point
GetSectionNumber returns the section number (for cyclic or periodic domains)
GetXmin returns the minimum of x-coordinates of the interpolation grid
GetXmax returns the maximum of x-coordinates of the interpolation grid
GetYmin returns the minimum of y-coordinates of the interpolation grid
GetYmax returns the maximum of y-coordinates of the interpolation grid
GetZmin returns the minimum of z-coordinates of the interpolation grid
GetZmax returns the maximum of z-coordinates of the interpolation grid
SetNbSubdivisions sets the number of points used for the intermediate regular grid
SetSubdivisionStep sets the spacings for the intermediate regular grid
SetNbCyclicSections sets the number of sections for a cyclic domain
SetNbCartesianSections sets the number of sections for a periodic domain
Append adds points that will be localized
InitInterpolationGrid initialisation of the intermediate regular grid for localization of points of a mesh
PreLocalizePoints first localization of points on the intermediate regular grid
LocalizePoints localization of points on the mesh

Public methods of GridInterpolationFull

GetOutputType returns the type of interpolation grid
SetOutputType sets the type of interpolation grid
GetNbPointsX returns the number of points along x
GetNbPointsY returns the number of points along y
GetNbPointsZ returns the number of points along z
SetNbPointsX sets the number of points along x
SetNbPointsY sets the number of points along y
SetNbPointsZ sets the number of points along z
GetCenter returns the intersection of the three planes
GetXmin returns the minimum of x-coordinates of the interpolation grid
GetXmax returns the maximum of x-coordinates of the interpolation grid
GetYmin returns the minimum of y-coordinates of the interpolation grid
GetYmax returns the maximum of y-coordinates of the interpolation grid
GetZmin returns the minimum of z-coordinates of the interpolation grid
GetZmax returns the maximum of z-coordinates of the interpolation grid
SetXmin sets the minimum of x-coordinates of the interpolation grid
SetXmax sets the maximum of x-coordinates of the interpolation grid
SetYmin sets the minimum of y-coordinates of the interpolation grid
SetYmax sets the maximum of y-coordinates of the interpolation grid
SetZmin sets the minimum of z-coordinates of the interpolation grid
SetZmax sets the maximum of z-coordinates of the interpolation grid
GetType returns the type of interpolation grid from a string
GetDimension returns the dimension of the output interpolation grid
GenerateGridPoints generates the points of the interpolation grid
InitGrid appends points to an object GridInterpolation in order to localize those points
Write writes header of the interpolation grid in binary format
WriteText writes header of the interpolation grid in text format
WriteVtk writes header of the interpolation grid in vtk format
Read reads header of the interpolation grid in binary format
ReadText reads header of the interpolation grid in text format
ReadVtk reads header of the interpolation grid in vtk format

Specific methods of GridInterpolationFull<Dimension2>

SetPlaneOutput regular points in 2-D
SetPlaneAxiOutput regular points in a 3-D plane
SetLineOutput regular points on a 2-D line
SetLineAxiOutput regular points on a 3-D line
SetPointOutput a single 2-D point
SetPointAxiOutput a single 3-D point
SetCircleOutput regular points on a 2-D circle
SetCircleAxiOutput regular points on a 3-D circle
SetPointsFileOutput 2-D points contained in a file
SetPointsFileAxiOutput 3-D points contained in a file
SetThreePlanesAxiOutput regular points on three cartesian planes in 3-D
SetVolumeAxiOutput regular points in 3-D

Specific methods of GridInterpolationFull<Dimension3>

SetPlaneOutput regular points in a 3-D plane
SetLineOutput regular points on a 3-D line
SetPointOutput a single 3-D point
SetCircleOutput regular points on a 3-D circle
SetPointsFileOutput 3-D points contained in a file
SetThreePlanesOutput regular points on three cartesian planes in 3-D
SetVolumeOutput regular points in 3-D

Methods of MeshInterpolation

GetNbSubdivisions returns the number of sub-intervals
SetRegularSubdivisions each edge of the mesh will be split into r regular subintervals
SetLobattoSubdivisions each edge of the mesh will be split into r subintervals, intermediate points are Gauss-Lobatto points
GetSubdivisionStep returns the intermediate points used for each edge
SetRegularPoints each edge of the mesh will be split into r regular subintervals
SetGaussQuadrature preparation of the computation of values of the solution on Gauss quadrature points
PointsReferenceVolume intermediate points on the unit tetrahedron/pyramid/prism and cube
PointsReferenceSurface intermediate points on the unit triangle/square
GetFileName returns the name of the file where the subdivided mesh will be stored
SetFileName sets the name of the file where the subdivided mesh will be stored
SetType sets the type of required output
GetType returns the type of required output
SetSurfaceReference sets the references of the surface for which solution will be interpolated
GetSurfaceReference returns the references of the surface for which solution will be interpolated
NumberingArray numbers for each element
InitProjectionVolume computes interpolation projector for volumetric interpolation
InitProjectionSurface computes interpolation projector for surfacic interpolation
ProjectVolume computes interpolation of the solution in an element
ProjectSurface computes interpolation of the solution on a surface of an element
IsVolumetric returns true if the object performs volumetric interpolations
IsSurfacic returns true if the object performs surfacic interpolations
ComputeSurfaceMesh extracts surface mesh and computes interpolation points
StoreGeometricDataSurface computes and stores interpolation points, normales, jacobian matrices, etc
GetElementNumberOfSurface returns element number of the element adjacent to the face
GetLocalPositionOfSurface returns local number of the face in the element
GetNbBoundary returns the number of faces on which interpolation will be computed
GetNbPointsQuadrature returns the number of quadrature points on all the surface mesh
PointsQuadrature returns the quadrature points of face i
Dfjm1Quadrature returns the jacobian matrices of face i
NormaleQuadrature returns the normales of face i
WeightsQuadrature returns the quadrature weights of face i
DsQuadrature returns the surface elements of integration of face i
GetNbAllQuadraturePoints returns the total number of quadrature points (on all the processors)
GetQuadraturePoint returns the global quadrature point (among all the processors)
GetQuadratureNormale returns the global normale (among all the processors)
GetQuadratureWeight returns the global quadrature weight (among all the processors)
ComputeEnHnOnBoundary computes u and du/dn (E x n and H x n for Maxwell's equastions) on quadratures points of the faces

Public attributes of VarOutputProblem

output_grid_param output parameters for interpolation on predefined grids
output_mesh_param output parameters for interpolation on subdivided meshes
all_points_display object for the localization of points of interpolation grids in the current mesh
var_grid object for interpolation on points of predefined grids
var_mesh_interp object for interpolation on subdivided mesh
grid_to_be_computed if true the points of interpolation grid need to be localized on the mesh

Public methods of VarOutputProblem

InitVarGrid initialisation of interpolation grids and outputs on subdivided meshes
ComputeVarGrid localization of points of predefined grids on the mesh
SetNbCyclicModes initialisation of the number of sections for output grids
InitOutput initialisation of outputs (before time iterations)
WriteSnapshot writes snapshots at a given time
ChangeTimeSnapshot changes the times associated for each output
InitVarMeshVolumetric computation of the subdivided mesh
InitVarMeshSurfacic computation of the subdivided mesh on a surface of the mesh
GetIndexOutputFiles returns the index in array output_grid_param for a type of predefined grid
WriteDatas writes outputs on predefined grids and subdivided meshes
GetEnHn forms u and du/dn (or E x n and H x n) from u and gradient of u
GetEnHnQuadrature forms u and du/dn (or E x n and H x n) from u and gradient of u
InitFftComputation initialisation of fft computations
GetModalOutput computes true value of the solution from its modal decomposition
ModifyOutputUnknown modification of the output vector u before writing it on the file
ModifyOutputGradient modification of the gradient vector du before writing it on the file
ComputeInterpolationU interpolation of u on points of predefined grids
ComputeInterpolationGradU interpolation of gradient of u on points of predefined grids
ComputeNodalUgradU interpolation of u and gradient of u on nodal points of the mesh
ProjectQuadratureToDof computation of u on degrees of freedom from values on quadrature points
WriteOutputFile writes the solution u on output files

Public methods of ParamOutputClass

GetComponent returns component of the solution to write
GetNature returns the type of interpolation grid associated with the output
SetNature sets the type of interpolation grid associated with the output
SetComponent sets component of the solution to write
IncrementSnapshot skips to next snapshot to write
GetSnapshotNumber returns the current snapshot number
GetTotalFieldFile returns the name of the file where the total field will be written
GetDiffractedFieldFile returns the name of the file where the diffracted field will be written
SetTotalFieldFile sets the name of the file where the total field will be written
SetDiffractedFieldFile sets the name of the file where the diffracted field will be written
GetTotalFieldFile returns the name of the file where the total field will be written
GetFileName returns the name of the i-th file
SetFileName sets the name of the i-th file
GetNbFile returns the name of extra-files handled by this object
InitTime initialisation of the snapshot number
SnapshotToStore returns true if a snapshot has to be written for this time
ChangeTime changes the current time
InitResult initialisation of output vectors
UpdateResult updating output vectors
WriteSnapshot writes snapshots on the disk
WriteNewSnapshot writes snapshots on the disk

affiche_mode

Syntax :

affiche_mode(root, extension, n1, n2)

This function displays all the eigenmodes whose number is comprised between n1 and n2. The extension .dat is automatically added, so extension may be empty.

affiche_mode('../totalSquareMode', '_U0', 0, 10);
affiche_mode('../totalSquareMode', '', 0, 10);

loadND, loadE, loadEH, loadH1

Syntax :

[X, Y, Z, coor, V] = loadND(fichier)
[X, Y, Z, coor, V] = loadND(fichier, num_inst, raff)
[X, Y, Z, coor, Ex, Ey, Ez] = loadE(base)
[X, Y, Z, coor, Ex, Ey, Ez, Hx, Hy, Hz] = loadEH(base)
[X, Y, Z, coor, U, dUx, dUy, dUz] = loadH1(base)

This function loads an output file of Montjoie. It only works if binary format (double or single precision) has been specified in Montjoie (this is actually the default format). num_inst is the snapshot number (equal to one if not specified), it is useful if you have specified outputs on differents planes for example (with two SismoPlane in the data file). raff is the level refinement (equal to one), it can be useful, if you wish to obtain a prettier solution (for example put 4 if you want to obtain a 400x400 solution from an output file containing 100x100 solution).

// loading a single field
[X, Y, Z, coor, V] = loadND('../totalCube_dU0.dat');

// if you want both u and its gradient (for a scalar field)
[X, Y, Z, coor, U, dUx, dUy, dUz] = loadE('../totalCube');

// if you are solving a vector field equation
// you can load simultaneously U0, U1 and U2
[X, Y, Z, coor, Ex, Ey, Ez] = loadE('../totalCube');

// for Maxwell's equation, it may be convenient
// to read both electric field E and magnetic field H
[X, Y, Z, coor, Ex, Ey, Ez, Hx, Hy, Hz] = loadEH('../totalCube');

// you can select another snapshot
[X, Y, Z, coor, V] = loadND('../totalCube_dU0.dat', 2);

// or refine the solution on more points
[X, Y, Z, coor, V] = loadND('../totalCube_dU0.dat', 1, 2);

contour2D

Syntax :

[non_loc, num_interface] = contour2D(X, Y, V);

This function returns the indexes for which the solution is null (it corresponds to points outside the domain), and the indexes at the boundary of the domain.

[X, Y, Z, coor, V] = loadND('../totalCarre_U0.dat');
[non_loc, num_inter ] = contour2D(X, Y, V);
// then you put nan on boundary, if you want to see the contour in white
V(num_inter) = nan;
plot2dinst(X, Y, real(V), -0.1, 0.1);

find_contour

Syntax :

[non_loc, num_interface] = find_contour(X, Y, Z,  V);

This function returns the indexes for which the solution is null (it corresponds to points outside the domain), and the indexes placed on the boundary of the domain.

[X, Y, Z, coor, V] = loadND('../totalCube_U0.dat');
[non_loc, num_inter ] = find_contour(X, Y, V);
// then you put nan on boundary, if you want to see the contour in
white
V(num_inter) = nan;
plot3d_planec(X, Y, Z, real(V), coor, -0.1, 0.1);

entier_to_string

Syntax :

s = entier_to_string(i)

This function converts an integer into a string such that it contains always four characters. If the integer is equal to 23, the corresponding string will be '0023'

s = entier_to_string(n);
[X, Y, Z, coor, V] = loadND(['../totalCarre', s, '_U0.dat']);

erreurL2

Syntax :

err = erreurL2(V, Vref);

This function returns the relative L2 error between two solutions. The two solutions are compared only on the points where they are both non-null (to avoid problems with the computation of interpolation grid).

[X, Y, Z, coor, V] = loadND('../totalCarreH02_U0.dat');
[X, Y, Z, coor, Vref] = loadND('../totalCarreH01_U0.dat');
err = erreurL2(V, Vref);

// you can compare as many components as you want
[X, Y, Z, coor, Ex, Ey, Ez] = loadE('../totalCarreH02');
[X, Y, Z, coor, ExRef, EyRef, EzRef] = loadE('../totalCarreH01');
err = erreurL2(Ex, ExRef, Ey, EyRef, Ez, EzRef);

erreurMediane

Syntax :

err = erreurMediane(V, Vref);

This function returns the relative "median" L2 error between two solutions. The two solutions are compared only on the points where they are both non-null (to avoid problems with the computation of interpolation grid). To avoid extrema values which are often due to problems of localization of the points in the mesh, this function actually computes the L2 by selecting the value for which 90% of points produce an error inferior to the error obtained with this value.

[X, Y, Z, coor, V] = loadND('../totalCarreH02_U0.dat');
[X, Y, Z, coor, Vref] = loadND('../totalCarreH01_U0.dat');
erreurMediane(V, Vref);

// you can compare as many components as you want
[X, Y, Z, coor, Ex, Ey, Ez] = loadE('../totalCarreH02');
[X, Y, Z, coor, ExRef, EyRef, EzRef] = loadE('../totalCarreH01');
erreurMediane(Ex, ExRef, Ey, EyRef, Ez, EzRef);

film2D

Syntax :

film2D(base, ext, debut, fin, cmin, cmax, draw_contour, num_inst, raff);
film2D(base, ext, debut, fin, cmin, cmax);

This function produces a picture .jpeg for each snapshot whose number is comprised between 'debut' and 'fin'. (cmin, cmax) defines the color scale. The optional argument draw_contour is equal to by default. If you put 1, the function contour2D will be used to display the boundary of the computational domain. The optional argument num_inst is equal to 1 by default. If you put another integer, it will read the corresponding snapshot contained in the file (it occurs if you have asked several SismoPlane).

// basic command
// it will read files '../TimeCarre0000_U1.dat', '../TimeCarre0001_U1.dat', etc
// and produces .jpeg for each snapshot
// jpeg files have the same root as output files, with .jpeg extension
film2D('../TimeCarre', '_U1', 0, 500, -0.001, 0.001);

// if you want to see the boundary
film2D('../TimeCarre', '_U1', 0, 500, -0.001, 0.001, 1);

film3D

Syntax :

film3D(base, ext, debut, fin, cmin, cmax, draw_contour, num_inst, raff);
film3D(base, ext, debut, fin, cmin, cmax);

This function produces a picture .jpeg for each snapshot whose number is comprised between 'debut' and 'fin'. (cmin, cmax) defines the color scale. The optional argument draw_contour is equal to by default. If you put 1, the function find_contour will be used to display the boundary of the computational domain. The optional argument num_inst is equal to 1 by default. If you put another integer, it will read the corresponding snapshot contained in the file (it occurs if you have asked several SismoPlane).

// basic command
// it will read files '../TimeCube0000_U1.dat', '../TimeCube0001_U1.dat', etc
// and produces .jpeg for each snapshot
// jpeg files have the same root as output files, with .jpeg extension
film3D('../TimeCube', '_U1', 0, 500, -0.001, 0.001);

// if you want to see the boundary
film3D('../TimeCube', '_U1', 0, 500, -0.001, 0.001, 1);

load1D, load1D_real, load1D_int

Syntax :

V = load1D(file_name)
V = load1D_real(file_name)
V = load1D_int(file_name)

These functions load a vector which has been previously written in a C++ code. For example, if you have written in C++, the following code :

Vector<double> V(10);
V.FillRand();
V.Write("sol_real.dat");

Vector<complex<double> > W(10);
W.FillRand();
W.Write("sol_complex.dat");

Vector<int> N(10);
N.FillRand();
N.Write("num.dat");

The following Python/Matlab code can be used :

V = load1D_real('sol_real.dat');
W = load1D('sol_complex.dat');
N = load1D_int('num.dat');

loadComplexMat

Syntax :

A = loadComplexMat(file_name)

This function loads a sparse complex matrix from a file. For example, if you have written in C++, the following code :

Matrix<complex<double>, General, ArrayRowSparse> A(10, 10);
A.FillRand();
A.Write("mat.dat");

The created file mat.dat will look like :

1 3 (0.5432143,-2.34433)
2 4 (2.90843,0.0)
3 2 (-1.232409,-0.090982)

The following Matlab code can be used to read the matrix :

A = loadComplexMat('mat.dat');

load_full, load_fullSym, load_fullReal, load_fullSym_real

Syntax :

A = load_full(file_name)
A = load_fullSym(file_name)
A = load_fullReal(file_name)
A = load_fullSym_real(file_name)

These functions load a dense matrix which has been previously written in a C++ code. For example, if you have written in C++, the following code :

Matrix<double> A(10, 10);
A.FillRand();
A.Write("Areal.dat");

Matrix<complex<double> > B(10, 10);
B.FillRand();
B.Write("Bcomplex.dat");

Matrix<double, Symmetric, RowSymPacked> C(10, 10);
C.FillRand();
C.Write("Csym.dat");

Matrix<complex<double>, Symmetric, RowSymPacked> D(10, 10);
D.FillRand();
D.Write("Dsym_complex.dat");

The following Python/Matlab code can be used to read those matrices :

A = load_fullReal('Areal.dat');
B = load_full('Bcomplex.dat');
C = load_fullSym_real('Csym.dat');
D = load_fullSym('Dsym_complex.dat');

plot_triangle, log_triangle

Syntax :

plot_triangle(position, order, taille, couleur)
log_triangle(position, order, taille, couleur)

These functions adds to the current figure a small triangle to check the order of convergence of a curve of the figure. You specify the position where the triangle will be placed, the expected order of convergence, the size of the triangle, and the color. log_triangle should be considered if you are using logarithmic scales (e.g. by using loglog command in Matlab).

// adding a magenta triangle at position [-1.2, -2.5] with a slope of 2
plot_triangle([-1.2, -2.5], 2, 0.2, 'm');

// for logarithmtic scales :
log_triangle([0.001, 1e-4], 2, 0.1, 'm');

plot2dinst

Syntax :

plot2dinst(X, Y, V);
plot2dinst(X, Y, V, cmin, cmax);

This function displays a 2-D solution. If the color extrema cmin, cmax are not specified, they are set to the minimum and maximum of the solution. With this function, you cannot plot a complex solution, you have to choose to display either its real part, imaginary or modulus.

// displaying imaginary part of a 2-D solution
[X, Y, Z, coor, V] = loadND('../totalCarre.dat');
plot2dinst(X, Y, imag(V));

plot3d_plane

Syntax :

plot3d_plane(X, Y, Z, V, coor);
plot3d_plane(X, Y, Z, V, coor, cmin, cmax);
plot3d_planec(X, Y, Z, V, coor, cmin, cmax);

This function displays a 3-D solution (on the three planes, whose intersection is coor). If the color extrema cmin, cmax are not specified, they are set to the minimum and maximum of the solution. With this function, you cannot plot a complex solution, you have to choose to display either its real part, imaginary or modulus. Use plot3d_planec to have the solution simultaneously on the three planes (in the top left corner).

// displaying modulus of a 3-D solution
[X, Y, Z, coor, V] = loadND('../totalCube.dat');
plot3d_planec(X, Y, Z, abs(V), coor);

str2num_python

Syntax :

str2num_python(epsilon)

This function converts a double into a string as it is performed by Python function str.

eps = 1.0;
s = str2num_python(eps);
[X, Y, Z, coor, V] = loadND(['total',s,'.dat']);

write_full_s, write_full_d

Syntax :

write_full_s(file_name)
write_full_d(file_name)

This function writes a dense matrix in binary format (it is equivalent to method Write for Seldon dense matrices). It is useful if you want to specify a variable index in Montjoie, you can create the index in Matlab, and use that function to write the data in a file, that Montjoie is able to read. For the name of the function, s means single precision, whereas d means double precision.

// creating a field
X = linspace(-5, 5, 200);
Y = linspace(-5, 5, 200);
[XI, YI] = meshgrid(X, Y);
V = exp(-(X.*X + Y.*Y))
// then writing V on the disk
write_full_d('rho.dat');

Medit

If in the data file, you have specified outputs on meshes (keywords SimoMeshSurfacic and SismoMeshVolumic), you can read the produced outputs with Medit. To display the solution, you have to type medit followed by the output file (extension.mesh) generated :

medit totalCarre_U0_real.mesh

Then type l (like Leopold) and m (like Mary) to see the solution. Look at the documentation of medit (you can type h to have the list of keyboard shortcuts displayed) for more details.