#include <stdio.h>
#include <stdlib.h>
#include <math.h>
#include <gsl/gsl_sf_bessel.h>
#include <gsl/gsl_integration.h>
#include <gsl/gsl_statistics_double.h>
#include "crusde_api.h"
Include dependency graph for pinel_thickplate.c:
Defines | |
#define | SIZE_W 100000 |
Functions | |
double | zeroth_order_bessel_integral (double e, void *params) |
double | first_order_bessel_integral (double e, void *params) |
double | f (double x, void *params) |
const char * | get_name () |
const char * | get_version () |
const char * | get_authors () |
const PluginCategory | get_category () |
const char * | get_description () |
void | register_parameter () |
Register parameters this Kernel claims from the input. | |
int | register_output_fields () |
Register output fields for spatial directions this Green's function calculates. | |
void | request_plugins () |
Request necessary plugins: crusde_request_green_plugin("elastic halfspace (pinel)"). | |
void | init () |
Initialization of the convolution. Allocation of memory for inputs and outputs. | |
void | run () |
Performs the fast convolution. | |
void | clear () |
Clean-up before this plug-in gets unloaded. | |
int | get_value_at (double **result, int x, int y, int t) |
Returns the Green's Function value at Point(x,y) at time 0. | |
Variables | |
double | g |
gravity [m/s^2] | |
double | E |
Young's Modulus [GPa]. | |
double | nu |
Poisson's ratio [-]. | |
double | rho_f |
Density of the inviscid fluid [kg/m^3]. | |
double | H |
thickness of the elastic plate [m] | |
double | U_vert_const |
Constant part of the vertical displacement equation. | |
double | U_hori_const |
Constant part of the horizontal displacement equation. | |
double | xx = 0.0 |
double | yy = 0.0 |
double | theta = 0.0 |
double | r = 0.0 |
double | local_r = 0.0 |
double | result_j0 = -1.0 |
double | result_j1 = -1.0 |
double | error = 0.0 |
int | quadrant = -1 |
short | sin_sign = -23 |
short | cos_sign = -23 |
double | one_minus_nu_squared = 0.0 |
double | H_squared = 0.0 |
double | H2 = 0.0 |
double | g_rho_f_E_term = 0.0 |
double | g_rho_f_E_term2 = 0.0 |
double | nu4 = 0.0 |
double | sin_h |
double | cos_h |
double | A |
double | B |
double | D |
double | global_e |
int | x_pos |
int | y_pos |
int | z_pos |
double | backpack [3] |
array for return values | |
double * | elastic = NULL |
keeps return values from elastic HS model | |
int(* | elastic_halfspace )(double **, int, int, int)=NULL |
gsl_integration_workspace * | w |
gsl_function | j0_integral |
gsl_function | j1_integral |
double | lower_integral_bound |
double | upper_integral_bound |
int | i = 0 |
Boundary conditions:<br /> .... The convention is such that a compression corresponds to a positive stress. At the lower surface of of the plate (at z=H) a pressure equal to the hydrostatic pressure in the fluid and no tangential component are applied. Displacements are calculated by integration of Bessel functions. At the surface (z=0) they are expressed by:
where
where c$ is the variable of integration. If
and
are Bessel functions of the first kind of zeroth- and first-order, respectively.
and
determine vertical and radial final relaxed displacement at point
, respectively. The coefficients
and
are defined by:
, then
and
which gives the solution of the elastic half-space.
The Hankel transform integral is calculated in the bounds using gsl_integration_qag of the GNU Scientific Library (http://www.gnu.org/software/gsl/).
|
|
|
Clean-up before this plug-in gets unloaded. Free the integration workspace! |
|
|
|
|
|
|
|
|
|
|
|
|
|
Returns the Green's Function value at Point(x,y) at time 0.
|
|
|
|
Initialization of the convolution. Allocation of memory for inputs and outputs. This function must not be called before register_parameter() unless none of the necessary values depends on parameters provided by the user, which are only set after they have been registered. This function is called some time after register_parameter(). Here U_vert_const and U_hori_const are initialized as follows:
|
|
Register output fields for spatial directions this Green's function calculates. This function calls crusde_register_output_field defined in crusde_api.h to register output field in the following order:
|
|
Register parameters this Kernel claims from the input. This function calls register_green_param() defined in crusde_api.h to register references to parameters this Green's function will need to operate properly. In case an XML is used to configure the experiment, the reference to this parameter will be identified by the string passed as second argument to register_green_param(). This function registers the references in the following order:
|
|
Request necessary plugins: crusde_request_green_plugin("elastic halfspace (pinel)"). empty |
|
Performs the fast convolution. writes data for actual time step to file |
|
|
|
|
|
|
|
array for return values
|
|
|
|
|
|
|
|
Young's Modulus [GPa].
|
|
keeps return values from elastic HS model
|
|
|
|
|
|
gravity [m/s^2]
|
|
|
|
|
|
|
|
thickness of the elastic plate [m]
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
Poisson's ratio [-].
|
|
|
|
|
|
|
|
|
|
|
|
|
|
Density of the inviscid fluid [kg/m^3].
|
|
|
|
|
|
|
|
Constant part of the horizontal displacement equation.
|
|
Constant part of the vertical displacement equation.
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|