Developers doc

Coupling elastic with other codes

Elastic is designed to be efficiently coupled to other simulation programs. For example, a meshing tool like mmg3d might regularize or adapt a mesh based on a displacement field computed by elastic.Or a finite element fluid solver might supply a nodal velocity field as boundary conditions,for elastic to extend this velocity field regularly into the domain.

Elastic can be used as a function in library, which can be called by another code. In this case, the other code is the master routine and calls elastic as needed. The function call has options that allow to be invoked with minimal overhead should you need to perform multiple runs within the same master program.

This section explains how to interface elastic in C/C++

Building elastic as a library

Elastic can be build as a dynamic (shared) library, which can be called from an application or a scripting language. It will be dynamically linked to the other program and loaded at run time. In a terminal, enter the directory LinearElasticity/build and type

make 

this command should produce a result like (on a Mac OS X):

Building C object CMakeFiles/Elas.dir/sources/elasti1_2d.c.o
[ 20%] Linking C shared library libElas.dylib
[ 80%] Built target Elas
[ 90%] Linking C executable elastic
[100%] Built target elastic

Then, to install the library and executable files, enter:

make install

and you get the following message (on a Mac OS X):

[ 80%] Built target Elas
[100%] Built target elastic
Install the project...
-- Install configuration: ""
-- Installing: /Users/logname/lib/libElas.dylib
-- Up-to-date: /Users/logname/include/elastic.h
-- Up-to-date: /Users/logname/include/ls_calls.h
-- Installing: /Users/logname/bin/elastic

(Note that for a shared library to be called by a program, all the auxiliary libraries it depends on must also exist as shared libraries, for instance libCommons.dylib or libCommons.so).

The library must be located at a place where the calling program can find it. By convention a shared library has a special name called soname. The soname has the prefix lib, the name of the library, the suffix .so, followed eventually by a period and a version number that is incremented whenever modifications are made. In Linux, the environnement variable

LD_LIBRARY_PATH

is a colon-separated set of directories where libraries should be searched for first, before the standard set of directories. Set this variable to the path where your elastic library has to be found.

Linking a program with elastic

You may be interested in called elastic as a library function in your own C/C++ code. This can be easily achieved by including the following header in your source code:

#include "ls_calls.h"

Here is an example of a main routine that calls elastic library functions:

int main(int argc,char *argv[]) {
  Mesh    mesh;
  LSst    *lsst;
  double  v[2];
  int     ier,np,na,nt,ne;

  memset(&mesh,0,sizeof(Mesh));                  /* default values */
  if ( !parsar(argc,argv,&mesh) )  return(1);    /* parse arguments in command line */

  fprintf(stdout,"  -- Test Elastic from library\n");
  fprintf(stdout,"\n - LOADING DATA\n");
  lsst = loadMesh(&mesh);                        /* load mesh file */
  if ( !lsst )  return(1);
  fprintf(stdout," - COMPLETED.\n");

  /* assign boundary conditions */
  fprintf(stdout,"\n ** MODULE TEST ELASTIC\n");
  LS_setPar(lsst,'+',0);                         /* verbosity set to '+' */
  v[0] = 0.0;   v[1] =  0.0;                     /* Dirichlet boundary conditions */
  LS_setBC(lsst,Dirichlet,1,'v',LS_ver,v);       /* on vertices labelled '1' */
  v[0] = 0.0;   v[1] = -0.2;
  LS_setBC(lsst,Dirichlet,2,'v',LS_ver,v);       /* on vertices labelled '2' */
  
  /* body force */
  v[0] = 0.0;   v[1] = -9.81;
  LS_setGra(lsst,v);                             /* set gravity vector */

  LS_setLame(lsst,0,186000.0,3400.0);            /* assign material properties */

  /* solve linear elasticity */
  if ( !LS_elastic(lsst) )  return(1);
  LS_stop(lsst);
  fprintf(stdout," ** COMPLETED\n\n");

  /* write out solution */
  fprintf(stdout," - WRITING DATA\n");
  if ( !saveSol(&mesh,lsst) )  return(1);
  fprintf(stdout," - COMPLETED\n");

  return(0);
}

The main program loads the mesh, assign boundary conditions, body forces and surface loads and call the module LS_elastic to solve the problem.

The function loadMesh must contain the following lines (in 2d):

LSst *loadMesh(pMesh mesh) {
  LSst     *lsst;
  double    c[2];
  int       k,i,int,ref,dim,ver,v[3];
  char      data[256];

  strcpy(data,mesh->name);
  if (!(inm = GmfOpenMesh(data,GmfRead,&ver,&dim)) ) {
    fprintf(stderr," # %s: file not found.\n",data);
    return(0);
  }
  fprintf(stdout,"    %s:",data);

  mesh->np = GmfStatKwd(inm,GmfVertices);
  mesh->nt = GmfStatKwd(inm,GmfTriangles);

  lsst = LS_init(dim,ver,P1,0);           /* memory allocation */ 
  LS_mesh(lsst,np,0,nt,0);    /* only vertices and triangles, no edges and tetra */

  /* read vertices */
  GmfGotoKwd(inm,GmfVertices);
  for (k=1; k<=np; k++) {
    GmfGetLin(inm,GmfVertices,&c[0],&c[1],&ref);  /* read coordinates and label */
    LS_addVer(lsst,k,c,ref);                      /* store vertex k in mesh */
  }
  /* read triangles */
  GmfGotoKwd(inm,GmfTriangles);
  for (k=1; k<=nt; k++) {
    GmfGetLin(inm,GmfTriangles,&v[0],&v[1],&v[2],&ref);
    LS_addTri(lsst,k,v,ref);
  }
  GmfCloseMesh(inm);
  fprintf(stdout," %d vertices, %d triangles\n",np,nt);

  return(lsst);
}

and the routine to save the solution

/* save solution */
int saveSol(pMesh mesh,LSst *lsst) {
  double      *u,buf[2];
  int          i,k,ia,outm,type,typtab[GmfMaxTyp];
  char        *ptr,data[128];

  strcpy(data,mesh->name);
  ptr = strstr(data,".mesh");
  if ( ptr )  {
    *ptr = '\0';
    strcat(data,".sol");
  }
  else {
    ptr = strstr(data,".sol");
    if ( !ptr )  strcat(data,".sol");
  }

  if ( !(outm = GmfOpenMesh(data,GmfWrite,mesh->ver,mesh->dim)) ) {
    fprintf(stderr," # unable to open %s\n",data);
    return(0);
  }
  fprintf(stdout,"    %s:",data);

  mesh->ver = GmfDouble;
  type = 1;
  typtab[0] = GmfVec;

  GmfSetKwd(outm,GmfSolAtVertices,mesh->np,type,typtab);
  u = LS_getSol(lsst);
  for (k=0; k<mesh->np; k++) {
    ia = mesh->dim*k;
    for (i=0; i<mesh->dim; i++)  buf[i] = u[ia+i];
    GmfSetLin(outm,GmfSolAtVertices,buf);
  }
  GmfCloseMesh(outm);
  fprintf(stdout," %d data vectors\n",mesh->np);

  return(1);
}