|
MPI_jacsolve
|
#include <vector>#include <list>#include <algorithm>#include <iostream>#include <fstream>#include <string>#include <cassert>#include "geom.h"Functions | |
| void | IniGeom (const int myid, const int procx, const int procy, int neigh[], int &color) |
| void | IniCoord (const int myid, const int procx, const int procy, float &xl, float &xr, float &yb, float &yt) |
| void | GetCoordsInRectangle (const int nx, const int ny, const float &xl, const float &xr, const float &yb, const float &yt, float xc[]) |
| void | GetConnectivityInRectangle (const int nx, const int ny, int ia[]) |
| void | GetBound (const int ib, const int nx, const int ny, const float w[], float s[]) |
| void | AddBound (const int ib, const int nx, const int ny, float w[], const float s[]) |
| void | SaveVectorP (const int myid, const string &name, const float u[], const int nx, const int ny, const float &xl, const float &xr, const float &yb, const float &yt) |
| void | GetMesh (const int nx, const int ny, const float xl, const float xr, const float yb, const float yt, int &nnode, float *&xc, int &nelem, int *&ia) |
| void | DebugMesh (const int nnode, const float xc[], const int nelem, const int ndof_e, const int ia[]) |
| void AddBound | ( | const int | ib, |
| const int | nx, | ||
| const int | ny, | ||
| float | w[], | ||
| const float | s[] | ||
| ) |
Computes w := w + s at the interface/boundary nodes on the boundary ib . South (ib==1), East (ib==2), North (ib==3), West (ib==4)
| [in] | ib | my local boundary |
| [in] | nx | number of discretization intervals in x-direction |
| [in] | ny | number of discretization intervals in y-direction |
| [in,out] | w | vector for all nodes of local discretization |
| [in] | s | short vector with values on boundary ib |
| void DebugMesh | ( | const int | nnode, |
| const float | xc[], | ||
| const int | nelem, | ||
| const int | ndof_e, | ||
| const int | ia[] | ||
| ) |
Prints the information for a finite element mesh
| [in] | nnode | number coordinates |
| [in] | xc[nnode] | array of coordinates |
| [in] | nelem | number of elements |
| [in] | ndof_e | degrees of freedom per element (= number of geometric points) |
| [in] | ia[nelem][ndof_e] | element connectivity array for linear triangular elements |
| void GetBound | ( | const int | ib, |
| const int | nx, | ||
| const int | ny, | ||
| const float | w[], | ||
| float | s[] | ||
| ) |
Copies the values of w corresponding to boundary ib onto vector s. South (ib==1), East (ib==2), North (ib==3), West (ib==4). The vector s has to be long enough!!
| [in] | ib | my local boundary |
| [in] | nx | number of discretization intervals in x-direction |
| [in] | ny | number of discretization intervals in y-direction |
| [in] | w | vector for all nodes of local discretization |
| [out] | s | short vector with values on boundary ib |
| void GetConnectivityInRectangle | ( | const int | nx, |
| const int | ny, | ||
| int | ia[] | ||
| ) |
Determines the element connectivity of linear triangular elements of a FEM discretization of a rectangle using nx x ny equidistant intervals for discretization.
| [in] | nx | number of discretization intervals in x-direction |
| [in] | ny | number of discretization intervals in y-direction |
| [out] | ia | element connectivity matrix with ia(3*s,3*s+1,3*s+2) as node numbers od element s |
| void GetCoordsInRectangle | ( | const int | nx, |
| const int | ny, | ||
| const float & | xl, | ||
| const float & | xr, | ||
| const float & | yb, | ||
| const float & | yt, | ||
| float | xc[] | ||
| ) |
Determines the coordinates of the dicretization nodes of the domain [xl, xr] x [yb, yt] which is discretized into nx x ny intervals.
| [in] | nx | number of discretization intervals in x-direction |
| [in] | ny | number of discretization intervals in y-direction |
| [in] | xl | x-coordinate of left boundary |
| [in] | xr | x-coordinate of right boundary |
| [in] | yb | y-coordinate of lower boundary |
| [in] | yt | y-coordinate of upper boundary |
| [out] | xc | coordinate vector of length 2n with x(2*k,2*k+1) as coodinates of node k |
| void GetMesh | ( | const int | nx, |
| const int | ny, | ||
| const float | xl, | ||
| const float | xr, | ||
| const float | yb, | ||
| const float | yt, | ||
| int & | nnode, | ||
| float *& | xc, | ||
| int & | nelem, | ||
| int *& | ia | ||
| ) |
The 2D domain [xl,xr] x [yb,yt] is discretized into a triangular mesh based on nx times ny equidistant grid of the rectangle.
The memory for coordinates array xc and for the element connectivity ia is allocated inside the function.
| [in] | nx | number of equidistant intervals in x-direction |
| [in] | ny | number of equidistant intervals in y-direction |
| [in] | xl | x-coordinate of left boundary |
| [in] | xr | x-coordinate of right boundary |
| [in] | yb | y-coordinate of lower boundary |
| [in] | yt | y-coordinate of upper boundary |
| [out] | nnode | number of coordinates |
| [out] | xc | array[nnode][2] of x- and y-coordinates |
| [out] | nelem | number of triangular elements |
| [out] | ia | array[nelem][3] element connectivity for linear triangular elements |
| void IniCoord | ( | const int | myid, |
| const int | procx, | ||
| const int | procy, | ||
| float & | xl, | ||
| float & | xr, | ||
| float & | yb, | ||
| float & | yt | ||
| ) |
The quadratic domain [0,1] x [0,1] is divided into procx * procy subdomains numbered rowise. According to my process with rank number myid, the coordinates of the lower left corner (xl, yb) and of the upper right corner (xr, yt) are generated.
| [in] | myid | my rank |
| [in] | procx | number of processes in x-direction |
| [in] | procy | number of processes in y-direction |
| [out] | xl | x-coordinate of left boundary |
| [out] | xr | x-coordinate of right boundary |
| [out] | yb | y-coordinate of lower boundary |
| [out] | yt | y-coordinate of upper boundary |
| void IniGeom | ( | const int | myid, |
| const int | procx, | ||
| const int | procy, | ||
| int | neigh[], | ||
| int & | color | ||
| ) |
The quadratic domain [0,1] x [0,1] is divided into procx * procy subdomains numbered rowise. According to my process with rank number myid, the vector neigh contains the ranks of the neighbours t the South, East, North, West.
| [in] | myid | my rank |
| [in] | procx | number of processes in x-direction |
| [in] | procy | number of processes in y-direction |
| [out] | neigh | vector of 4 elements containing the ranks of neighbouring processes, -1 indicates a boundary in this direction |
| [out] | color | red/black coloring of subdomains (0/1) |
| void SaveVectorP | ( | const int | myid, |
| const string & | name, | ||
| const float | u[], | ||
| const int | nx, | ||
| const int | ny, | ||
| const float & | xl, | ||
| const float & | xr, | ||
| const float & | yb, | ||
| const float & | yt | ||
| ) |
Stores the values of vector u of subdomain myid into a file. The file stores rowise the x- and y- coordinates together with the value from u . The domain [xl, xr] x [yb, yt] is discretized into nx x ny intervals.
| [in] | myid | my rank |
| [in] | name | basename of file name (file name will be extended by the rank number) |
| [in] | u | local vector |
| [in] | nx | number of discretization intervals in x-direction |
| [in] | ny | number of discretization intervals in y-direction |
| [in] | xl | x-coordinate of left boundary |
| [in] | xr | x-coordinate of right boundary |
| [in] | yb | y-coordinate of lower boundary |
| [in] | yt | y-coordinate of upper boundary |
1.8.11