Difference between revisions of "User Defined Functions"

Latest revision as of 15:56, 6 November 2023

Some commands let you set properties using a function that is evaluated at run time.

Function Variables

When a user-defined function option is allowed in any command, you can enter any valid function of the following variables when doing MPM simulations:

x - x position in length units
y - y position in length units
z - z position in length units
t - current time in alt time units
dt - the time step in alt time units
q - particle rotation in radians about the z (or θ if axisymmetric) axis

For FEA calculations, the following variables are allowed and may refer to position of element centroid or node depending on the command:

x - x position in length units
y - y position in length units
r - radial position of axisymmetric calculations in length units (synonym for x).
z - axial position of axisymmetric calculations in length units (synonym for x)

When a user-defined function is used, it will be calculated using these variables and should return a result in the units expected by the command. Note that commands that allow functions may only allow a subset of these variables (due to command context). You can refer to each command for the allowed variables. For example, some MPM options (as detailed in their documentation) require the function to depend only on time. Particle-based MPM boundary conditions let the function depend on clockwise particle rotational angle q about the z axis (in radians), which allows rotation of the boundary conditions with the particle. Note that q is particle rotation since the start of the simulation and will differ from material angle if the particle started with a non-zero orientation angle (i.e., the current material angle is the sum of q and its initial angle).

When setting up MPM simulations or FEA calculations, a few more variables are sometimes allowed:

R - radial position in length units for axisymmetric calculation; R is a synonym for x, which also works
Z - axial position in length units for axisymmetric calculation; Z is a synonym for y, which also works
D - distance from origin length units
T - counter-clockwise angle (in radians) from the positive x axis (i.e., θ in polar coordinates)

Here (R,Z) are axisymmetric coordinates and (D,T) are polar coordinates, where D is distance from the origin to the (x,y) (or (R,Z) if axisymmetric) point (in length units) and T is counter-clockwise angle (in radians) from the positive x (or R if axisymmetric) axis to the point. These extra variables are only allowed in:

MPM Region commands or <Body> blocks to set initial particle velocity, angular momentum, or rotation angle (also for development membranes in OSParticulas.
<Deform> block to set initial particle displacement.
FEA Area command or <Area> block to set initial material angle.
FEA Temperature command or <Temperature> element to set initial temperature on the nodes.

Function Format

Some details on entering functions are:

In scripted files, the function must be enclosed in quotes (e.g., "x^2+y^2") to prevent it from being evaluated as a command expression prior to being used in the analysis.
Enter variables simply as x, y, etc., and not with the preceding "#" used for command expression variables.
Operators: The function uses standard operators + - * / and ^ with standard operator precedence for addition, subtraction, multiplication, division, and raising to a power.
The function can contain the following defined functions:
- sin(x), cos(x), and tan(x) - trigonometric function of angle in radians.
- asin(x), acos(x), and atan(x) - inverse trigonometric function with result in radians.
- sqrt(x) - square root.
- log(x) and ln(x) - log base 10 and natural log, respectively (note these are different than log functions used in command expressions in scripted files due to different math expression parser in the code engines).
- exp(x) - exponential of x.
- abs(x) - absolute value of x.
- int(x) - integer part of x as next lower integer. For negative numbers, it is next lower integer or int(-3.4)=-4.
- Sinh(x), Cosh(x), and Tanh(x) - hyperbolic trigonometric functions. (note the initial uppercase is needed here, and differs from scripted command expressions due to different math expression parser in the code engines).
- erf(x) and erfc(x) - error function and error function complement.
The following functions are helpful if creating boundary conditions:
- [math]\displaystyle{ {\rm ramp}(A,x) = \left\{\begin{array}{ll} 0 & x\le 0 \\ Ax & 0\lt x\lt 1 \\ A & x \ge 1 \end{array}\right. }[/math]
- [math]\displaystyle{ {\rm cosramp}(A,x) = \left\{\begin{array}{ll} 0 & x\le 0 \\ {A\over 2}\left(1-\cos(\pi x)\right) & 0\lt x\lt 1 \\ A & x \ge 1 \end{array}\right. }[/math]
- [math]\displaystyle{ {\rm box}(A,x) = \left\{\begin{array}{ll} 0 & x\le 0 \\ A & 0\lt x\lt 1 \\ 0 & x \ge 1 \end{array}\right. }[/math]
- [math]\displaystyle{ {\rm sinbox}(A,x) = \left\{\begin{array}{ll} 0 & x\le 0 \\ A\sin(\pi x) & 0\lt x\lt 1 \\ 0 & x \ge 1 \end{array}\right. }[/math]
- [math]\displaystyle{ {\rm sgn}(x) = \left\{\begin{array}{ll} -1 & x \lt 0 \\ 0 & x = 0 \\ +1 & x \gt 0 \end{array}\right. }[/math]
- [math]\displaystyle{ {\rm sign}(x) = \left\{\begin{array}{ll} 0 & x \lt 0 \\ 1 & x \ge0 \end{array}\right. }[/math]
- [math]\displaystyle{ {\rm tri}(x) = \left\{\begin{array}{ll} 1-|x| & |x| \lt 1 \\ 0 & |x| \ge1 \end{array}\right. }[/math]
- [math]\displaystyle{ {\rm mod}(x,y) = {\rm fmod}(x,y) = x-{\rm floor}(x/y)*y }[/math]
For example, to apply a ramp from 0 to t_f, use x = t/t_f; to apply a ramp from t_s to t_f,use x = (t-t_s)/(t_f-t_s).
Functions can include pi (or Pi or PI) for the number π.
Exponential Notation: numbers can have "e" or "E" for powers of ten such as 1.4e3 for 1400.
Extra spaces in the function are ignored.

Notes

Integrals and Derivatives

The cosramp(A,x) and sinbox(A,x) are useful for gradually applying velocity or displacement at the start of the simulation rather than instantaneously applying some velocity or displacement at time zero. To interchange between velocity and displacement, it is useful to evaluate their integrals and derivatives:

[math]\displaystyle{ \int_0^x {\rm sinbox}(A,u)du = {A\over\pi}\bigl(1-\cos(\pi x)\bigr) = {\rm cosramp}\left({2A\over\pi},x\right) }[/math]

[math]\displaystyle{ \int_0^x {\rm cosramp}(A,u)du = \left\{\begin{array}{ll} {A\over 2}\left(x-{\sin(\pi x)\over \pi}\right) & x\lt 1 \\ A\left(x-{1\over 2}\right) & x \ge 1 \end{array}\right. = {A\over 2}\bigl(x + (x-1){\rm sign}(x-1)\bigr) - {\rm sinbox}\left({A\over2\pi},x\right) }[/math]

[math]\displaystyle{ {d\bigl({\rm cosramp}(A,x)\bigr)\over dx} = {A\pi\over 2}\sin(\pi x) = {\rm sinbox}\left({A\pi\over 2},x\right) }[/math]

[math]\displaystyle{ {d\bigl({\rm sinbox}(A,x)\bigr)\over dx} = \left\{\begin{array}{ll} 0 & x\le0 \\ {A\pi}\cos(\pi x) & 0\lt x\lt 1 \\ 0 & x \ge 1 \end{array}\right. = \pi A\thinspace {\rm sign}(x) - {\rm cosramp}(2\pi A,x) }[/math]

Cosine Ramp to Desired Velocity

In problems with constant velocity loading, it might be preferable to ramp velocity to the desired velocity V over some ramp time t_r rather start the problem at time zero with velocity V. The ramped-up velocity can be done using the cosramp(A,x) function by setting A = V and x = t/t_r or:

[math]\displaystyle{ v(t) = {\rm cosramp}\left(V,{t\over t_r}\right) }[/math]

Integrating this velocity, the resulting displacement is:

[math]\displaystyle{ d(t) = \left\{\begin{array}{ll} {V\over 2}\left(t-{t_r\over \pi}\sin\left({\pi t\over t_r}\right)\right) & t\lt t_r \\ V\left(t-{t_r\over 2}\right) & t \ge t_r \end{array}\right. }[/math]

In terms of defined functions, this displacement can be written

[math]\displaystyle{ d(t) = {V\over2}\left[t + (t-t_r){\rm sign}(t-t_r) - {\rm sinbox}\left({t_r\over\pi},{t\over t_r}\right)\right] }[/math]

For example, to ramp to V and then hold constant velocity until reaching desired displacement d_max, the simulation time should be set to:

[math]\displaystyle{ t_{max} = {d_{max}\over V} + {t_r\over 2} }[/math]

The first term is the simulation time using constant velocity while the second term is the extra time needed because the velocity was ramped to velocity V over time t_r.

Cosine Ramp to Desired Displacement

Rather that apply displacement instantaneously, which would be an impact problem, it might be preferable to ramp to displacement d over ramp time t_r. A displacement function would be

[math]\displaystyle{ d(t) = {\rm cosramp}\left(d,{t\over t_r}\right) }[/math]

Differentiating this displacement to get velocity that must be applied as a boundary condition gives

[math]\displaystyle{ v(t) = {\pi d\over 2t_r}\sin\left({\pi t\over t_r}\right) = {\rm sinbox}\left({\pi d\over 2t_r},{t\over t_r}\right) }[/math]

The maximum velocity used during this loading is

[math]\displaystyle{ v_{max} = {\pi d\over 2t_r} }[/math]

Linear Ramp to Desired Displacement

To ramp to a displacement d over time t_r with a constant velocity, the ramp velocity is

[math]\displaystyle{ v(t) = {\rm box}\left({d\over t_r},{t\over t_r}\right) }[/math]

The displacement corresponding to this velocity is

[math]\displaystyle{ d(t) = {\rm ramp}\left(d,{t\over t_r}\right) }[/math]

Loading and Unloading Cycles

To ramp at velocity v_r from t₀ and t₀+t_r then unload using the reverse velocity from t₀+t_r to t₀+2t_r, the velocity is:

[math]\displaystyle{ v(t) = {\rm box}\left(v_r,{t-t_0\over t_r}\right)-{\rm box}\left(v_r,{t-t_0-t_r\over t_r}\right) = {\rm box}\left(v_r,{t-t_0\over 2t_r}\right){\rm sgn}(t_0+t_r-t) }[/math]

The displacement corresponding to this velocity is

[math]\displaystyle{ d(t) = v_rt_r {\rm tri}\left({t-t_0\over t_r}-1\right) }[/math]

A sequence of loading and unloading steps can be created by defining periodic function as explained in the next section.

Periodic Boundary Condition Function

The mod(A,x) function is useful for defining any periodic boundary condition function of time. First create any function of time that is well defined over a given cycle time t_c. To define a periodic function replace every "t" in that function with "mod(t,t_c)". The result will be a periodic function that will repeat the same defined cycle as long as the simulations keeps running. For example, a sequence of triangular displacements from the previous section's example staring at the beginning (t₀=0) could be defined with:

[math]\displaystyle{ v(t) = {\rm box}\left(v_r,{{\rm mod}(t,2t_r)\over t_r}\right)-{\rm box}\left(v_r,{{\rm mod}(t,2t_r)-t_r\over t_r}\right) }[/math]

Rigid Rotation

See Rigid Material documentation for details on moving a plot of rigid particles in rigid rotation about a point.