|
| 1 | +/*! |
| 2 | + * \file IntegrateTest4-c.c |
| 3 | + * \brief |
| 4 | + * \author Thomas Helfer |
| 5 | + * \date 21/09/2018 |
| 6 | + * \copyright (C) Copyright Thomas Helfer 2018. |
| 7 | + * Use, modification and distribution are subject |
| 8 | + * to one of the following licences: |
| 9 | + * - GNU Lesser General Public License (LGPL), Version 3.0. (See accompanying |
| 10 | + * file LGPL-3.0.txt) |
| 11 | + * - CECILL-C, Version 1.0 (See accompanying files |
| 12 | + * CeCILL-C_V1-en.txt and CeCILL-C_V1-fr.txt). |
| 13 | + */ |
| 14 | + |
| 15 | +#include <math.h> |
| 16 | +#include <stdio.h> |
| 17 | +#include <stdlib.h> |
| 18 | +#include "MGIS/ThreadPool.h" |
| 19 | +#include "MGIS/Model/Model.h" |
| 20 | +#include "MGIS/Behaviour/MaterialDataManager.h" |
| 21 | +#include "MGIS/Behaviour/Integrate.h" |
| 22 | + |
| 23 | +int test_status = EXIT_SUCCESS; |
| 24 | +mgis_model_Model* model = NULL; |
| 25 | +mgis_ThreadPool* p = NULL; |
| 26 | +mgis_bv_MaterialDataManager* m = NULL; |
| 27 | + |
| 28 | +static void check_status(const mgis_status s) { |
| 29 | + if (s.exit_status != MGIS_SUCCESS) { |
| 30 | + fprintf(stderr, "invalid function call: %s\n", s.msg); |
| 31 | + mgis_model_free_model(&model); |
| 32 | + mgis_bv_free_material_data_manager(&m); |
| 33 | + mgis_free_thread_pool(&p); |
| 34 | + exit(EXIT_FAILURE); |
| 35 | + } |
| 36 | +} // end of check_status |
| 37 | + |
| 38 | +static int check(const int b, const char* const e) { |
| 39 | + if (b == 0) { |
| 40 | + test_status = EXIT_FAILURE; |
| 41 | + fprintf(stderr, "%s\n", e); |
| 42 | + } |
| 43 | + return b; |
| 44 | +} // end of check |
| 45 | + |
| 46 | +int main(const int argc, const char* const* argv) { |
| 47 | + if (check(argc == 2, "expected two arguments") == 0) { |
| 48 | + return EXIT_FAILURE; |
| 49 | + } |
| 50 | + const mgis_real eps = 1e-10; |
| 51 | + const mgis_real dt =0.1; |
| 52 | + mgis_size_type o; |
| 53 | + mgis_real* isvs0; /* internal state variables at the beginning of the time step */ |
| 54 | + mgis_real* isvs1; /* internal state variables at the end of the time step */ |
| 55 | + mgis_size_type isvs_stride; /* internal state variables stride */ |
| 56 | + const mgis_size_type n = 100; |
| 57 | + mgis_bv_MaterialStateManager* |
| 58 | + s0; /* state at the beginning of the time step */ |
| 59 | + mgis_bv_MaterialStateManager* s1; /* state at the end of the time step */ |
| 60 | + mgis_real xi[11]; /* values of 'x' for the first integration point */ |
| 61 | + mgis_real xe[11]; /* values of 'x' for the last integration point */ |
| 62 | + mgis_size_type ni,ne; |
| 63 | + mgis_size_type idx; |
| 64 | + mgis_size_type i; |
| 65 | + mgis_real A; |
| 66 | + mgis_real t; |
| 67 | + mgis_real x_ref; |
| 68 | + int r; |
| 69 | + check_status(mgis_create_thread_pool(&p, 2)); |
| 70 | + check_status(mgis_model_load(&model, argv[1], "ode_rk54", "Tridimensional")); |
| 71 | + check_status(mgis_bv_behaviour_get_parameter_default_value(&A, model, "A")); |
| 72 | + check_status(mgis_bv_create_material_data_manager(&m, model, 100)); |
| 73 | + check_status( |
| 74 | + mgis_bv_behaviour_get_internal_state_variable_offset(&o, model, "x")); |
| 75 | + check_status(mgis_bv_material_data_manager_get_state_0(&s0, m)); |
| 76 | + check_status(mgis_bv_material_data_manager_get_state_1(&s1, m)); |
| 77 | + /* initialize the internal state variable */ |
| 78 | + check_status( |
| 79 | + mgis_bv_material_state_manager_get_internal_state_variables(&isvs0, s0)); |
| 80 | + check_status( |
| 81 | + mgis_bv_material_state_manager_get_internal_state_variables(&isvs1, s1)); |
| 82 | + check_status( |
| 83 | + mgis_bv_material_state_manager_get_internal_state_variables_stride( |
| 84 | + &isvs_stride, s1)); |
| 85 | + for (idx = 0; idx != n; ++idx) { |
| 86 | + isvs1[idx * isvs_stride + o] = 1; |
| 87 | + } |
| 88 | + /* initialize the external state variable */ |
| 89 | + check_status( |
| 90 | + mgis_bv_material_state_manager_set_uniform_external_state_variable( |
| 91 | + s1, "Temperature", 293.15)); |
| 92 | + /* copy s1 in s0 */ |
| 93 | + check_status(mgis_bv_update_material_data_manager(m)); |
| 94 | + // integration */ |
| 95 | + ni = o; |
| 96 | + ne = (n - 1) * isvs_stride + o; |
| 97 | + xi[0] = isvs0[ni]; |
| 98 | + xe[0] = isvs0[ne]; |
| 99 | + for (i = 0; i != 10; ++i) { |
| 100 | + check_status(mgis_bv_integrate_material_data_manager( |
| 101 | + &r, p, m, MGIS_BV_INTEGRATION_NO_TANGENT_OPERATOR, dt)); |
| 102 | + check_status(mgis_bv_update_material_data_manager(m)); |
| 103 | + xi[i + 1] = isvs1[ni]; |
| 104 | + xe[i + 1] = isvs1[ne]; |
| 105 | + } |
| 106 | + check_status(mgis_model_free_model(&model)); |
| 107 | + check_status(mgis_bv_free_material_data_manager(&m)); |
| 108 | + check_status(mgis_free_thread_pool(&p)); |
| 109 | + t = 0; |
| 110 | + for (i = 0; i != 11; ++i) { |
| 111 | + x_ref = exp(-A * t); |
| 112 | + if (fabs(xi[i] - x_ref) > eps) { |
| 113 | + fprintf(stderr, |
| 114 | + "IntegrateTest: invalid value for x " |
| 115 | + "at the first integration point" |
| 116 | + "(expected '%g', computed '%g')\n", |
| 117 | + x_ref, xi[i]); |
| 118 | + return EXIT_FAILURE; |
| 119 | + } |
| 120 | + if (fabs(xe[i] - x_ref) > eps) { |
| 121 | + fprintf(stderr, |
| 122 | + "IntegrateTest: invalid value for x " |
| 123 | + "at the first integration point" |
| 124 | + "(expected '%g', computed '%g')\n", |
| 125 | + x_ref, xe[i]); |
| 126 | + return EXIT_FAILURE; |
| 127 | + } |
| 128 | + t += dt; |
| 129 | + } |
| 130 | + return EXIT_SUCCESS; |
| 131 | +} |
0 commit comments