User defined extensions
For many simulations the functions implemented in numgeo are sufficient. However, for more complex problems it may be necessary to extend the existing functionality. For this purpose, numgeo offers the possibility to connect own developments to predefined interfaces and to extend the functionality of numgeo. The following interfaces are available:
-
Amplitude: This subroutine is used for user-defined amplitude values, characterized by the parameter
*Amplitude, type=user
command in the input file. The header and variable description are given in the following. -
Boundary conditions: This subroutine is used for user-defined boundary values, characterized by the parameter:
*UBoundary
command in the input file. -
Initial stress: This subroutine is used for user-defined initial stresses, characterized by the parameter
user
on the*initial conditions, type = stress, user
command in the input file. The subroutine is called at the start of the analysis for each applicable material calculation point in the model and can be used to define all active initial stress components at material points as functions of coordinates, element number, integration point number, etc. -
Material model: This subroutine is used to extend numgeo with user-defined constitutive material models, characterized by
*Mechanical = user
command in the input file. The header and variable description is given in the following. -
Initial state variables: This subroutine is used to initialize user-defined state variables (
statev
) at particular material points (integration points), characterized by the parameteruser
on the*initial conditions
command in the input file:*initial conditions, type = state variables, user
-
Contact properties: This subroutine is used for user-defined contact properties, characterized by the parameter
user
on the*interaction, user
command in the input file. The subroutine is called in each calculation step and has to define the contact properties of the chosen constitutive contact model. -
Loads: This subroutine is used for user-defined distributed loads, characterized by the parameter
*Dload, user
command in the input file. The header and variable description are given in the following. The coordinates of all nodes of the underlying element are provided via the interface. However, only the nodes at the element face being load are assigned a value, even if a value is given to the other nodes of the element.
The following steps are recommended if you wish to use user routines for any of the above functions:
- Copy the corresponding sample user file
.f90
and corresponding bash file (.sh
on linux and.bat
on Windows) shipped with the numgeo executable to the directory where you want to run the simulation. - Modify the sample file to your needs, make sure to follow Fortran conventions, declare all newly introduced variables appropriately and never ever change the interface.
- Specify the use of the required user routine in the input file of your numgeo simulation
- Start the simulation
Precompiled user routines
- If you have a already compiled shared object of your user routine
.so
in linux and.dll
in Windows you can also place this file in the directory instead of the source file and bash file. - If the source file and the bash file of a user defined routine are present in the directory, numgeo will always recompile the user routine and override existing shared objects.
- It is always good to start from an input file that does not rely on user routines. If this file runs without troubles, incrementally add the required user routines.
user routines on Windows
To link numgeo to user subroutines on Windows (using Intel OneAPI), the installation of Intel Visual Studio is required. See Installation manual
Python Script for plotting the results
The following python script can be used to plot the results of the different simulations using user-routines as well as the basic simulation of the Rigid strip footing.
#%%
import numpy as np
import matplotlib.pyplot as plt
plt.rcParams['xtick.direction'] = 'in'
plt.rcParams['ytick.direction'] = 'in'
plt.rcParams['text.usetex'] = True
plt.rcParams['font.size'] = 12
plt.rcParams['axes.spines.right'] = False
plt.rcParams['axes.spines.top'] = False
def cm2inch(*tupl):
inch = 2.54
if isinstance(tupl[0], tuple):
return tuple(i/inch for i in tupl[0])
else:
return tuple(i/inch for i in tupl)
node_set_name = "FOUNDATION.TOP_LEFT_NODE_node_7703.dat"
ref = np.genfromtxt("./reference/simulation-print-out/"+node_set_name, skip_header=1)
#
# Reference solution
#
plt.figure(figsize=cm2inch(12,8))
plt.plot(ref[:,0], ref[:,2], label='linear load')
plt.xlabel('Time in s')
plt.ylabel('Settlement in m')
plt.legend(loc='lower left', frameon=False)
plt.tight_layout()
plt.savefig('./time_settlement_reference.pdf', dpi=400)
plt.savefig('./time_settlement_reference.tiff', dpi=400)
plt.savefig('./time_settlement_reference.png', dpi=400)
#
# User amplitude
#
node_set_name = "FOUNDATION.TOP_LEFT_NODE_node_7703"
data = np.genfromtxt("./user-amplitude-mesh-new/simulation-print-out/"+node_set_name, skip_header=1)
plt.figure(figsize=cm2inch(12,8))
plt.plot(ref[:,0], ref[:,2], label='linear load')
plt.plot(data[:,0], data[:,2], label='$t (sin(\omega_1 t)/5 - sin(\omega_2 t)/20)$')
plt.xlabel('Time in s')
plt.ylabel('Settlement in m')
plt.legend(loc='lower left', frameon=False)
plt.tight_layout()
plt.savefig('./time_settlement_user_amplitude.pdf', dpi=400)
plt.savefig('./time_settlement_user_amplitude.tiff', dpi=400)
plt.savefig('./time_settlement_user_amplitude.png', dpi=400)
#
# User material
#
data = np.genfromtxt("./user material/simulation-print-out/"+node_set_name, skip_header=1)
plt.figure(figsize=cm2inch(12,8))
plt.plot(ref[:,0], ref[:,2], label='$E=10$ MPa')
plt.plot(data[:,0], data[:,2], label='$E(p) = E_0 \\left( \\frac{-\sigma_{ii}}{3p_{ref}} \\right)^n$ \n with: $E_0=5$ MPa, $n=0.2$, $p_{ref}=1$ kPa')
plt.xlabel('Time in s')
plt.ylabel('Settlement in m')
plt.legend(loc='lower left', frameon=False)
plt.tight_layout()
plt.savefig('./time_settlement_user_material.pdf', dpi=400)
plt.savefig('./time_settlement_user_material.tiff', dpi=400)
plt.savefig('./time_settlement_user_material.png', dpi=400)
data = np.genfromtxt("./user material/simulation-print-out/"+node_set_name, skip_header=1)
plt.figure(figsize=cm2inch(12,8))
plt.plot(data[:,0], data[:,2], label='$E(p) = E_0 \\left( \\frac{-\sigma_{ii}}{3p_{ref}} \\right)^n$ \n with: $E_0=5$ MPa, $n=0.2$, $p_{ref}=1$ kPa')
plt.xlabel('Time in s')
plt.ylabel('Settlement in m')
plt.legend(loc='lower left', frameon=False)
plt.tight_layout()
plt.savefig('./time_settlement_user_material2.pdf', dpi=400)
plt.savefig('./time_settlement_user_material2.tiff', dpi=400)
plt.savefig('./time_settlement_user_material2.png', dpi=400)
#
# User initial state
#
data = np.genfromtxt("./user material/simulation-print-out/"+node_set_name, skip_header=1)
plt.figure(figsize=cm2inch(12,8))
# plt.plot(ref[:,0], ref[:,2], label='$E=10$ MPa')
plt.plot(data[:,0], data[:,2], label='$E(p) = E_0 f_p$')
data = np.genfromtxt("./user initial state/simulation-print-out/"+node_set_name, skip_header=1)
plt.plot(data[:,0], data[:,2], label='$E(p) = E_0 f_p f_e$ with $f_e=\\frac{1+e}{e}$ \n ...and linear decrease of $e_{0}(d)$')
plt.xlabel('Time in s')
plt.ylabel('Settlement in m')
plt.legend(loc='lower left', frameon=False)
plt.tight_layout()
plt.savefig('./time_settlement_user_initial_state.pdf', dpi=400)
plt.savefig('./time_settlement_user_initial_state.tiff', dpi=400)
plt.savefig('./time_settlement_user_initial_state.png', dpi=400)
#
# User unsat material
#
data = np.genfromtxt("./user material unsat/simulation-print-out/"+node_set_name, skip_header=1)
plt.figure(figsize=cm2inch(12,8))
plt.semilogx(data[:,0], data[:,2], label='$E(p) = E_0 \\left( 1 + \\alpha s (S^e)^\\beta \\right)$ \n with: $E_0=7.5$ MPa, $\\alpha=1.0$, $\\beta=1.0$')
plt.xlabel('Time in s')
plt.ylabel('Settlement in m')
plt.legend(loc='upper right', frameon=False)
plt.tight_layout()
plt.savefig('./time_settlement_user_material_unsat.pdf', dpi=400)
plt.savefig('./time_settlement_user_material_unsat.tiff', dpi=400)
plt.savefig('./time_settlement_user_material_unsat.png', dpi=400)