This header defines specific integrators functions in C. More...

#include "mbs_data.h"
#include "mbs_aux.h"

Typedefs
typedef void(*	initialize_integrator_ptr) (MbsData mbs_data, MbsDirdyn mbs_dd)
	These three pointers of functions should be defined for each integrator separately, in the corresponding files They are called in dirdyn, according to the integrator chosen. More...

typedef int(*	loop_integrator_ptr) (double t0, double tf, MbsData mbs_data, MbsDirdyn mbs_dd)
	main loop of integration, in which the integrator is called More...

typedef void(*	finish_integrator_ptr) (MbsData mbs_data, MbsDirdyn mbs_dd)
	end of the integration, free the memory More...

Enumerations
enum	{ RK4, Dopri5, Rosenbrock, EulerEx, Eulaire, EulerIm, Bader, WMethods, AlphaM, Custom }
	Enumeration of the possible integrator for dirdyn. More...

Functions
void	set_integrator (MbsDirdyn *dd)
	Set the function pointer in MbsDirdyn. More...

void	print_warnings_integrator (MbsData mbs_data, MbsDirdyn mbs_dd, int type_of_integrator)
	Check the options set by user to warn him when he modified an unused option for the integrator. More...

void	print_warnings_constant_step_integrator (MbsDirdyn mbs_dd, char integrator_name)
	Print the warning message, called by print_warnings_integrator for constant step size dt integrators. More...

int	rk4 (double y[], double dydx[], int n, double x, double h, double yout[], int(derivs)(double, double[], double[], MbsData , MbsDirdyn ), MbsData s, MbsDirdyn *dd)
	Runge Kutta 4 integrator implementation Given values for the variables y[1..n] and their derivatives dydx[1..n] known at x, use the fourth-order Runge-Kutta method to advance the solution over an interval h and return the incremented variables as yout[1..n], which need not be a distinct array from y. More...

int	rosenbrock (int n, int(derivs)(double, double[], double[], MbsData , MbsDirdyn ), double x, double y[], double eps, double hmax, double htry, long nmax, double dydx[], double yscal[], double hnext, MbsData s, MbsDirdyn dd, double hdid)
	Fourth-order Rosenbrock step for integrating stiff problems, with monitoring of local truncation error to adjust stepsize. More...

int	ThetaSC (double y[], double dydx[], int n, double x, double h, double yout[], int(derivs)(double, double[], double[], MbsData , MbsDirdyn ), MbsData s, MbsDirdyn *dd)
	Unknown integrator. More...

int	mbs_estim_jac_acc (double x, double htry, double y[], double dydx[], int compute_dfdx, double dfdx[], double *dfdy, int n, int(derivs)(double, double[], double[], MbsData , MbsDirdyn ), MbsData s, MbsDirdyn dd)
	Evaluate the Jacobian of a function accelerations using finite difference. More...

int	mbs_freeze_jac (int freeze_index, int next_freeze_index, double x, double h, double y[], double dydx[], int compute_dfdx, double dfdx[], double dfdy, double dydx_freeze[], double dfdx_freeze[], double dfdy_freeze, int n, int(derivs)(double, double[], double[], MbsData , MbsDirdyn ), MbsData s, MbsDirdyn mbs_dd)
	Freeze the Jacobian of an integrator structure until mbs_dd->options->n_freeze. More...

int	euler_implicit (double y[], double dydx[], int n, double x, double h, double yout[], int(derivs)(double, double[], double[], MbsData , MbsDirdyn ), MbsData s, MbsDirdyn mbs_dd)
	Euler Implicit integrator implementation. More...

int	w_methods (double y[], double dydx[], int n, double x, double h, double yout[], int(derivs)(double, double[], double[], MbsData , MbsDirdyn ), MbsData s, MbsDirdyn *mbs_dd)
	W Methods integrator implementation. More...

Detailed Description

This header defines specific integrators functions in C.

Creation date: May 2018

Author: Sebastien Timmermans, Olivier Lantsoght

(c) Universite catholique de Louvain

Typedef Documentation

◆ finish_integrator_ptr

typedef void(* finish_integrator_ptr) (MbsData *mbs_data, MbsDirdyn *mbs_dd)

end of the integration, free the memory

◆ initialize_integrator_ptr

typedef void(* initialize_integrator_ptr) (MbsData *mbs_data, MbsDirdyn *mbs_dd)

These three pointers of functions should be defined for each integrator separately, in the corresponding files They are called in dirdyn, according to the integrator chosen.

initialization of the integrator structure(s)

◆ loop_integrator_ptr

typedef int(* loop_integrator_ptr) (double t0, double tf, MbsData *mbs_data, MbsDirdyn *mbs_dd)

main loop of integration, in which the integrator is called

Enumeration Type Documentation

◆ anonymous enum

anonymous enum

Enumeration of the possible integrator for dirdyn.

// Set the integrator in main : mbs_dirdyn->options->integrator = Dopri5; mbs_dirdyn->options->resfilename = "dirdyn_Dopri5";

Adding a custom integrator can be done using the "Custom" keyword and modifying the files mbs_custom.c / .h

Enumerator
RK4
Dopri5
Rosenbrock
EulerEx
Eulaire
EulerIm
Bader
WMethods
AlphaM
Custom

Function Documentation

◆ euler_implicit()

int euler_implicit	(	double	y[],
		double	dydx[],
		int	n,
		double *	x,
		double	h,
		double	yout[],
		int()(double, double[], double[], MbsData , MbsDirdyn *)	derivs,
		MbsData *	s,
		MbsDirdyn *	mbs_dd
	)

Euler Implicit integrator implementation.

@source (adapted from) Arnold M. et al, Linearly implicit time integration methods in real-time applications: DAEs and stiff ODEs Multibody System Dynamics, 2007, 17:99â117

Parameters

y	state vector of variable
dydx	derivative of y by x
n	number of variables
x	current time value, at which the function needs to be computed
h	time step
yout	solution of the incremented variables
derivs	the function computing f'
s	the MbsData structure of the model on which dirdyn analysis is computed.
mbs_dd	the MbsDirdyn structure related to the integration.

Returns: error_code Return a negative number if something went wrong, 0 otherwise

calls an external function

A = (- h* Jv - h*h * Jp )

A = ( I - * Jv - h * h * Jp)

add 1/h in diagonal terms (identity matrix)

Compute the product : Jp * Xv

B = qdd + h * Jp * Xv

delta p = Xv + h * delta v

with delta v = B

Xv (n+1) with delta v = B

Xp (n+1)

◆ mbs_estim_jac_acc()

int mbs_estim_jac_acc	(	double	x,
		double	htry,
		double	y[],
		double	dydx[],
		int	compute_dfdx,
		double	dfdx[],
		double **	dfdy,
		int	n,
		int()(double, double[], double[], MbsData , MbsDirdyn *)	derivs,
		MbsData *	s,
		MbsDirdyn *	dd
	)

Evaluate the Jacobian of a function accelerations using finite difference.

Modify and do not restore dydx. WARNING: All vector/matrix indexes start at 0! Update the value of dfdx and dfdy.

Parameters

x	current time value, at which the function needs to be computed
htry	time step
y	state vector of variable
compute_dfdx	flag to ON to compute the dfdx terms, default = 1 ; if compute_dfdx = 0 => dfdx is filled with zeroes
dydx	derivative of y by x
dfdx	derivative of f by x
dfdy	derivative of f by y
n	number of variables
derivs	the function computing f'
s	the MbsData structure of the model on which dirdyn analysis is computed.
dd	the MbsDirdyn structure related to the integration.

Returns: Error status, <0 in case of failure.

◆ mbs_freeze_jac()

int mbs_freeze_jac	(	int	freeze_index,
		int *	next_freeze_index,
		double	x,
		double	h,
		double	y[],
		double	dydx[],
		int	compute_dfdx,
		double	dfdx[],
		double **	dfdy,
		double	dydx_freeze[],
		double	dfdx_freeze[],
		double **	dfdy_freeze,
		int	n,
		int()(double, double[], double[], MbsData , MbsDirdyn *)	derivs,
		MbsData *	s,
		MbsDirdyn *	mbs_dd
	)

Freeze the Jacobian of an integrator structure until mbs_dd->options->n_freeze.

WARNING: All vector/matrix indexes start at 0!

Parameters

freeze_index	current count of refreezing steps (if 0, recompute Jac)
next_freeze_index	[out] pointor to next refreezing step
x	current time value, at which the function needs to be computed
h	time step
y	state vector of variable
dydx	derivative of y by x
compute_dfdx	flag to ON to compute the dfdx terms, default = 1 ; if compute_dfdx = 0 => dfdx is filled with zeroes
dfdx	derivative of f by x
dfdy	derivative of f by y
dydx_freeze	pointor to freezed derivative of y by x
dfdx_freeze	pointor to freezed derivative of f by x
dfdy_freeze	pointor to freezed derivative of f by y
n	number of variables
derivs	the function computing f'
s	the MbsData structure of the model on which dirdyn analysis is computed.
mbs_dd	the MbsDirdyn structure related to the integration.

Returns: Error status, <0 in case of failure.

◆ print_warnings_constant_step_integrator()

void print_warnings_constant_step_integrator	(	MbsDirdyn *	mbs_dd,
		char *	integrator_name
	)

Print the warning message, called by print_warnings_integrator for constant step size dt integrators.

Parameters

mbs_dd	MbsDirdyn structure
integrator_name	string to print in the warning message, corresponding to the integrator (RK4, EulerEx, ...)

◆ print_warnings_integrator()

void print_warnings_integrator	(	MbsData *	mbs_data,
		MbsDirdyn *	mbs_dd,
		int	type_of_integrator
	)

Check the options set by user to warn him when he modified an unused option for the integrator.

Parameters

mbs_data	MbsData structure
mbs_dd	MbsDirdyn structure
type_of_integrator	ID if the used integrator (RK4, Dopri5...)

◆ rk4()

int rk4	(	double	y[],
		double	dydx[],
		int	n,
		double	x,
		double	h,
		double	yout[],
		int()(double, double[], double[], MbsData , MbsDirdyn *)	derivs,
		MbsData *	s,
		MbsDirdyn *	dd
	)

Runge Kutta 4 integrator implementation Given values for the variables y[1..n] and their derivatives dydx[1..n] known at x, use the fourth-order Runge-Kutta method to advance the solution over an interval h and return the incremented variables as yout[1..n], which need not be a distinct array from y.

The user supplies the routine derivs(x,y,dydx), which returns derivatives dydx at x.

Parameters

y	state vector of variable
dydx	derivative of y by x
n	number of variables
x	current time value, at which the function needs to be computed
h	time step
yout	solution of the incremented variables
derivs	The function computing f'
s	the MbsData structure of the model on which dirdyn analysis is computed.
dd	the MbsDirdyn structure related to the integration.

Returns: Error status, <0 in case of failure.

◆ rosenbrock()

int rosenbrock	(	int	n,
		int()(double, double[], double[], MbsData , MbsDirdyn *)	derivs,
		double *	x,
		double	y[],
		double	eps,
		double	hmax,
		double	htry,
		long	nmax,
		double	dydx[],
		double	yscal[],
		double *	hnext,
		MbsData *	s,
		MbsDirdyn *	dd,
		double *	hdid
	)

Fourth-order Rosenbrock step for integrating stiff problems, with monitoring of local truncation error to adjust stepsize.

Input are the dependent variable vector y[1..n] and its derivative dydx[1..n] at the starting value of the independent variable x. WARNING: this feature has not been heavily tested !

@source: H. H. Rosenbrock, "Some general implicit processes for the numerical solution of differential equations", The Computer Journal (1963) 5(4): 329-330 Shampine, L.F. 1982, ACM Transactions on Mathematical Software, vol. 8, pp. 93â113

Parameters

n	Dimension of y
derivs	The function computing f'
x	Value of x
y	Initial y values, index starting at 0
eps	The required accuracy
hmax	The maximum stepsize allowed
htry	The stepsize to be attempted
nmax	The maximal number of allowed steps to reach tolerances
dydx	Initial y derivative values, index starting at 0
yscal	The vector [0..n-1] against which the error is scaled, index starting at 0
hnext	The estimated next stepsize
s	The MbsData structure of the model on which dirdyn analysis is computed.
dd	The MbsDirdyn structure related to the integration.
hdid	The stepsize that was actually accomplished

Returns: error_code Return a negative number if something went wrong, 0 otherwise

calls an external function

◆ set_integrator()

void set_integrator ( MbsDirdyn * dd )

Set the function pointer in MbsDirdyn.

Parameters

dd	MbsDirdyn structure to adapt

◆ ThetaSC()

int ThetaSC	(	double	y[],
		double	dydx[],
		int	n,
		double	x,
		double	h,
		double	yout[],
		int()(double, double[], double[], MbsData , MbsDirdyn *)	derivs,
		MbsData *	s,
		MbsDirdyn *	dd
	)

Unknown integrator.

Returns: Error status, <0 in case of failure.

◆ w_methods()

int w_methods	(	double	y[],
		double	dydx[],
		int	n,
		double	x,
		double	h,
		double	yout[],
		int()(double, double[], double[], MbsData , MbsDirdyn *)	derivs,
		MbsData *	s,
		MbsDirdyn *	mbs_dd
	)

W Methods integrator implementation.

@source (adapted from) Arnold M. et al, Linearly implicit time integration methods in real-time applications: DAEs and stiff ODEs Multibody System Dynamics, 2007, 17:99â117

Parameters

y	state vector of variable
dydx	derivative of y by x
n	number of variables
x	current time value, at which the function needs to be computed
h	time step
yout	solution of the incremented variables
derivs	the function computing f'
s	the MbsData structure of the model on which dirdyn analysis is computed.
mbs_dd	the MbsDirdyn structure related to the integration.

Returns: error_code Return a negative number if something went wrong, 0 otherwise

calls an external function

A = (- h* Jv - h*h * Jp )

A = ( I - * Jv - h * h * Jp)

add 1/h in diagonal terms (identity matrix)

partial

useless if s<2

sum alpha + gamma

Jp * Xv

Complete right-hand term B

B = f + h * Jp * Xv gamma +

delta p = Xv + h * delta v

Typedefs

Enumerations

Functions

Detailed Description

Typedef Documentation

◆ finish_integrator_ptr

◆ initialize_integrator_ptr

◆ loop_integrator_ptr

Enumeration Type Documentation

◆ anonymous enum

Function Documentation

◆ euler_implicit()

◆ mbs_estim_jac_acc()

◆ mbs_freeze_jac()

◆ print_warnings_constant_step_integrator()

◆ print_warnings_integrator()

◆ rk4()

◆ rosenbrock()

◆ set_integrator()

◆ ThetaSC()

◆ w_methods()