FEDEM Solver  R8.0
Source code of the dynamics solver
solverInterface.h
Go to the documentation of this file.
1 /* SPDX-FileCopyrightText: 2023 SAP SE
2  *
3  * SPDX-License-Identifier: Apache-2.0
4  *
5  * This file is part of FEDEM - https://openfedem.org
6  */
21 #ifndef _SOLVER_INTERFACE_H
22 #define _SOLVER_INTERFACE_H
23 
24 #include <stddef.h>
25 
26 #ifdef __cplusplus
27 extern "C" {
28 #else
29 #include <stdbool.h>
30 #endif
31 
65  int solverInit(int argc, char** argv, const char* cfsi = NULL,
66  const double* stateData = NULL, const int ndat = 0,
67  const double* gagesData = NULL, const int ngda = 0,
68  const double* extInput = NULL, int* nxinp = NULL);
69 
87  int restartFromState(const double* stateData, const int ndat,
88  const int writeToRDB = 2);
89 
110  bool solveWindow(int nStep, int nInc, int nIn, int nOut, const int* fOut,
111  const double* times, const double* inputs, double* outputs,
112  int nDat, double* stateData, int* ierr);
113 
117  int haveResults();
118 
124 
130 
137 
144 
150 
161  bool saveState(double* stateData, const int ndat);
162 
168  bool saveTransformationState(double* stateData, const int ndat);
169 
175  bool saveGages(double* stateData, const int ndat);
176 
185  bool getStrainsFromDisp(const double* disp, const int* gageIDs, double* eps,
186  const int ndof, const int ng);
187 
196  bool getBeamForcesFromDisp(const double* disp, const int* beamIDs,
197  double* forces, const int ndof, const int nBeam);
198 
207  bool getRelDisp(const double* disp, const int* Ids, double* relDis,
208  const int ndof, const int nId);
209 
218  bool getRespVars(const double* disp, const int* Ids, double* resVar,
219  const int ndof, const int nId);
220 
231  bool savePartDeformationState(int bid, double* stateData, const int ndat);
232 
243  bool savePartStressState(int bid, double* stateData, const int ndat);
244 
252  bool solveNext(int* ierr);
253 
267  bool startStep(int* ierr);
268 
283  bool solveIteration(int* ierr, bool all = false);
284 
298  bool solveEigenModes(const int nmod, double* eval, double* evec,
299  bool dofOrder, int useLaPack, int* ierr);
300 
314  bool solveInverse(const double* x,
315  const int* xeqs, const int* feqs,
316  int ndis, int nfrs, int* ierr);
317 
328  int solverDone(bool removeSingletons = true);
329 
335  void solverClose();
336 
352  int setExtFunc(int funcId, double value);
353 
363  double getTime(int tFlag, int* ierr);
364 
369  bool setTime(double nextTime);
370 
374  inline double getCurrentTime() { int ierr; return getTime(0,&ierr); }
375 
382  inline double getNextTime(int* ierr) { return getTime(1,ierr); }
383 
396  double evalFunc(int uid, const char* tag = NULL,
397  double x = -1.0, int* ierr = NULL);
398 
406  double getFunc(int uid, int* ie = NULL) { return evalFunc(uid,NULL,-1.0,ie); }
407 
411  int getFuncId(const char* tag);
412 
420  int getEquations(int bid, int* meqn);
421 
433  int getStateVar(int bid, double* var);
434 
439  int getSystemSize(bool dofs = false);
440 
446  bool getSystemMatrix(double* Smat, int iMat = 0);
447 
452  bool getNewtonMatrix(double* Nmat);
453 
458  bool getStiffnessMatrix(double* Kmat);
459 
464  bool getMassMatrix(double* Mmat);
465 
470  bool getDampingMatrix(double* Cmat);
471 
478  bool getElementStiffnessMatrix(double* Kmat, int bid);
479 
485  bool getRhsVector(double* Rvec, int iVec = 0);
486 
491  bool setRhsVector(const double* Rvec);
492 
497  bool addRhsVector(const double* Rvec);
498 
504  bool getJointSprCoeff(double* sprCoeff, int bid);
505 
512  const char* getFileName(const char* fileOpt, char* fileName, int nc);
513 
514 #ifdef __cplusplus
515 }
516 #endif
517 
518 #endif
integer, save ndof
Definition: diffractionModule.f90:16
integer(ptr), save, private x
Definition: extCtrlSysRoutinesModule.f90:16
double getNextTime(int *ierr)
Returns the physical time of the next step of the simulation.
Definition: solverInterface.h:382
bool getNewtonMatrix(double *Nmat)
Returns current content of the system Newton matrix.
bool savePartDeformationState(int bid, double *stateData, const int ndat)
Saves the deformation of an FE part to the provided core array.
int getStateVar(int bid, double *var)
Returns the current state variables for a specified object.
int getSystemSize(bool dofs=false)
Returns the dimension (number of equations) of the system.
int getTransformationStateSize()
Returns the required length of the state vector which contains translation and rotation information f...
bool solveIteration(int *ierr, bool all=false)
Solves current linearized equation system and updates the state.
bool getElementStiffnessMatrix(double *Kmat, int bid)
Returns the content of an element stiffness matrix.
bool solveInverse(const double *x, const int *xeqs, const int *feqs, int ndis, int nfrs, int *ierr)
Solves the inverse problem at current time/load step.
int haveResults()
Returns whether current time step have results to be saved.
bool getStiffnessMatrix(double *Kmat)
Returns current content of the system stiffness matrix.
bool getRhsVector(double *Rvec, int iVec=0)
Returns current content of the system right-hand-side vector.
bool getStrainsFromDisp(const double *disp, const int *gageIDs, double *eps, const int ndof, const int ng)
Gets gage strain values from given displacement field.
bool addRhsVector(const double *Rvec)
Updates the content of the system right-hand-side vector.
bool getBeamForcesFromDisp(const double *disp, const int *beamIDs, double *forces, const int ndof, const int nBeam)
Gets internal sectional forces from given displacements.
int restartFromState(const double *stateData, const int ndat, const int writeToRDB=2)
Restarts a simulation process from the provided state.
int solverInit(int argc, char **argv, const char *cfsi=NULL, const double *stateData=NULL, const int ndat=0, const double *gagesData=NULL, const int ngda=0, const double *extInput=NULL, int *nxinp=NULL)
Initializes a new solver process.
bool getSystemMatrix(double *Smat, int iMat=0)
Returns current content of the system Newton matrix.
double evalFunc(int uid, const char *tag=NULL, double x=-1.0, int *ierr=NULL)
Evaluates a specified general function in the model.
void solverClose()
Cleans up singleton objects in core memory.
bool savePartStressState(int bid, double *stateData, const int ndat)
Saves the stress state of an FE part to the provided core array.
int setExtFunc(int funcId, double value)
Assigns a sensor value from a physical twin to the simulation model.
double getFunc(int uid, int *ie=NULL)
Returns the current value of the specified general function.
Definition: solverInterface.h:406
int getFuncId(const char *tag)
Returns the user ID of a tagged function.
int solverDone(bool removeSingletons=true)
Closes down current model and cleans up core memory and on disk.
bool getDampingMatrix(double *Cmat)
Returns current content of the system damping matrix.
bool setTime(double nextTime)
Sets the time increment size of the next time step.
int getPartStressStateSize(int bid)
Returns the required length of the state vector which contains von Mises stress information for the s...
int getEquations(int bid, int *meqn)
Returns the equation numbers for the DOFs of a specified object.
bool startStep(int *ierr)
Starts a new time/load step.
bool getMassMatrix(double *Mmat)
Returns current content of the system mass matrix.
int getStateSize()
Returns the required length of the state vector to be is used when restarting a simulation from an in...
bool getRelDisp(const double *disp, const int *Ids, double *relDis, const int ndof, const int nId)
Gets relative distances between triads from given displacements.
int getGagesSize()
Returns the required length of the initial strain gage vector which is used when restarting a simulat...
double getTime(int tFlag, int *ierr)
Returns the physical time of a time step.
bool saveTransformationState(double *stateData, const int ndat)
Saves current transformation state to the provided core array.
bool solveNext(int *ierr)
Advances the solution one time/load step forward.
bool saveGages(double *stateData, const int ndat)
Saves the initial gage strains to the provided core array.
bool setRhsVector(const double *Rvec)
Replaces current content of the system right-hand-side vector.
int getPartDeformationStateSize(int bid)
Returns the required length of the state vector which contains deformation information for the specif...
bool solveWindow(int nStep, int nInc, int nIn, int nOut, const int *fOut, const double *times, const double *inputs, double *outputs, int nDat, double *stateData, int *ierr)
Executes the simulation for a specified time window.
bool saveState(double *stateData, const int ndat)
Saves current state to the provided array.
const char * getFileName(const char *fileOpt, char *fileName, int nc)
Gets a file name option from the command-line arguments.
bool getRespVars(const double *disp, const int *Ids, double *resVar, const int ndof, const int nId)
Gets response variable values from given displacements.
double getCurrentTime()
Returns the current physical time of the simulation.
Definition: solverInterface.h:374
bool solveEigenModes(const int nmod, double *eval, double *evec, bool dofOrder, int useLaPack, int *ierr)
Solves the eigenvalue problem at current state.
bool getJointSprCoeff(double *sprCoeff, int bid)
Gets spring stiffness coefficients for a joint.