Commit 9cf3737b authored by Kévin Moreau's avatar Kévin Moreau

Last commit using the "old" fast marching interface.

- Cleanup unused symbolic link
- Add "mode" atomic test case
- Pass pull_out test case from ring geometry to plain disk geometry
- Correct problems in infinite_dam_notch test case
parent 746ec9cd
Pipeline #112 skipped
......@@ -13,7 +13,7 @@ file(GLOB LIST
# levelset_update
# mesh_generation
# mesh_generation_2
# mode
mode
# negative
# nonlocal_projection
# patch_test
......
enable_testing()
add_test(hansbo_enrichment ${CMAKE_CURRENT_BINARY_DIR}/hansbo_enrichment ${CMAKE_CURRENT_BINARY_DIR}/data/square.msh)
Point(1) = {0, 0, 0, 1};
Point(2) = {1, 0, 0, 1};
Point(3) = {1, 1, 0, 1};
Point(4) = {0, 1, 0, 1};
Line(1) = {1, 2};
Line(2) = {2, 3};
Line(3) = {3, 4};
Line(4) = {4, 1};
Line Loop(5) = {3, 4, 1, 2};
Plane Surface(6) = {5};
Transfinite Line {1, 2, 3, 4} = 11 Using Progression 1;
This diff is collapsed.
$NOD
121
1 0.000000 0.000000 0.
2 0.100000 0.000000 0.
3 0.200000 0.000000 0.
4 0.300000 0.000000 0.
5 0.400000 0.000000 0.
6 0.500000 0.000000 0.
7 0.600000 0.000000 0.
8 0.700000 0.000000 0.
9 0.800000 0.000000 0.
10 0.900000 0.000000 0.
11 1.000000 0.000000 0.
12 0.000000 0.100000 0.
13 0.100000 0.100000 0.
14 0.200000 0.100000 0.
15 0.300000 0.100000 0.
16 0.400000 0.100000 0.
17 0.500000 0.100000 0.
18 0.600000 0.100000 0.
19 0.700000 0.100000 0.
20 0.800000 0.100000 0.
21 0.900000 0.100000 0.
22 1.000000 0.100000 0.
23 0.000000 0.200000 0.
24 0.100000 0.200000 0.
25 0.200000 0.200000 0.
26 0.300000 0.200000 0.
27 0.400000 0.200000 0.
28 0.500000 0.200000 0.
29 0.600000 0.200000 0.
30 0.700000 0.200000 0.
31 0.800000 0.200000 0.
32 0.900000 0.200000 0.
33 1.000000 0.200000 0.
34 0.000000 0.300000 0.
35 0.100000 0.300000 0.
36 0.200000 0.300000 0.
37 0.300000 0.300000 0.
38 0.400000 0.300000 0.
39 0.500000 0.300000 0.
40 0.600000 0.300000 0.
41 0.700000 0.300000 0.
42 0.800000 0.300000 0.
43 0.900000 0.300000 0.
44 1.000000 0.300000 0.
45 0.000000 0.400000 0.
46 0.100000 0.400000 0.
47 0.200000 0.400000 0.
48 0.300000 0.400000 0.
49 0.400000 0.400000 0.
50 0.500000 0.400000 0.
51 0.600000 0.400000 0.
52 0.700000 0.400000 0.
53 0.800000 0.400000 0.
54 0.900000 0.400000 0.
55 1.000000 0.400000 0.
56 0.000000 0.500000 0.
57 0.100000 0.500000 0.
58 0.200000 0.500000 0.
59 0.300000 0.500000 0.
60 0.400000 0.500000 0.
61 0.500000 0.500000 0.
62 0.600000 0.500000 0.
63 0.700000 0.500000 0.
64 0.800000 0.500000 0.
65 0.900000 0.500000 0.
66 1.000000 0.500000 0.
67 0.000000 0.600000 0.
68 0.100000 0.600000 0.
69 0.200000 0.600000 0.
70 0.300000 0.600000 0.
71 0.400000 0.600000 0.
72 0.500000 0.600000 0.
73 0.600000 0.600000 0.
74 0.700000 0.600000 0.
75 0.800000 0.600000 0.
76 0.900000 0.600000 0.
77 1.000000 0.600000 0.
78 0.000000 0.700000 0.
79 0.100000 0.700000 0.
80 0.200000 0.700000 0.
81 0.300000 0.700000 0.
82 0.400000 0.700000 0.
83 0.500000 0.700000 0.
84 0.600000 0.700000 0.
85 0.700000 0.700000 0.
86 0.800000 0.700000 0.
87 0.900000 0.700000 0.
88 1.000000 0.700000 0.
89 0.000000 0.800000 0.
90 0.100000 0.800000 0.
91 0.200000 0.800000 0.
92 0.300000 0.800000 0.
93 0.400000 0.800000 0.
94 0.500000 0.800000 0.
95 0.600000 0.800000 0.
96 0.700000 0.800000 0.
97 0.800000 0.800000 0.
98 0.900000 0.800000 0.
99 1.000000 0.800000 0.
100 0.000000 0.900000 0.
101 0.100000 0.900000 0.
102 0.200000 0.900000 0.
103 0.300000 0.900000 0.
104 0.400000 0.900000 0.
105 0.500000 0.900000 0.
106 0.600000 0.900000 0.
107 0.700000 0.900000 0.
108 0.800000 0.900000 0.
109 0.900000 0.900000 0.
110 1.000000 0.900000 0.
111 0.000000 1.000000 0.
112 0.100000 1.000000 0.
113 0.200000 1.000000 0.
114 0.300000 1.000000 0.
115 0.400000 1.000000 0.
116 0.500000 1.000000 0.
117 0.600000 1.000000 0.
118 0.700000 1.000000 0.
119 0.800000 1.000000 0.
120 0.900000 1.000000 0.
121 1.000000 1.000000 0.
$ENDNOD
$ELM
240
1 1 11 1 2 1 2
2 1 11 1 2 2 3
3 1 11 1 2 3 4
4 1 11 1 2 4 5
5 1 11 1 2 5 6
6 1 11 1 2 6 7
7 1 11 1 2 7 8
8 1 11 1 2 8 9
9 1 11 1 2 9 10
10 1 11 1 2 10 11
11 1 12 2 2 11 22
12 1 12 2 2 22 33
13 1 12 2 2 33 44
14 1 12 2 2 44 55
15 1 12 2 2 55 66
16 1 12 2 2 66 77
17 1 12 2 2 77 88
18 1 12 2 2 88 99
19 1 12 2 2 99 110
20 1 12 2 2 110 121
21 1 13 3 2 111 112
22 1 13 3 2 112 113
23 1 13 3 2 113 114
24 1 13 3 2 114 115
25 1 13 3 2 115 116
26 1 13 3 2 116 117
27 1 13 3 2 117 118
28 1 13 3 2 118 119
29 1 13 3 2 119 120
30 1 13 3 2 120 121
31 1 14 4 2 1 12
32 1 14 4 2 12 23
33 1 14 4 2 23 34
34 1 14 4 2 34 45
35 1 14 4 2 45 56
36 1 14 4 2 56 67
37 1 14 4 2 67 78
38 1 14 4 2 78 89
39 1 14 4 2 89 100
40 1 14 4 2 100 111
41 2 101 6 3 1 2 13
42 2 101 6 3 1 13 12
43 2 101 6 3 2 3 13
44 2 101 6 3 3 14 13
45 2 101 6 3 3 4 15
46 2 101 6 3 3 15 14
47 2 101 6 3 4 5 15
48 2 101 6 3 5 16 15
49 2 101 6 3 5 6 17
50 2 101 6 3 5 17 16
51 2 101 6 3 6 7 17
52 2 101 6 3 7 18 17
53 2 101 6 3 7 8 19
54 2 101 6 3 7 19 18
55 2 101 6 3 8 9 19
56 2 101 6 3 9 20 19
57 2 101 6 3 9 10 21
58 2 101 6 3 9 21 20
59 2 101 6 3 10 11 21
60 2 101 6 3 11 22 21
61 2 101 6 3 12 13 23
62 2 101 6 3 13 24 23
63 2 101 6 3 13 14 25
64 2 101 6 3 13 25 24
65 2 101 6 3 14 15 25
66 2 101 6 3 15 26 25
67 2 101 6 3 15 16 27
68 2 101 6 3 15 27 26
69 2 101 6 3 16 17 27
70 2 101 6 3 17 28 27
71 2 101 6 3 17 18 29
72 2 101 6 3 17 29 28
73 2 101 6 3 18 19 29
74 2 101 6 3 19 30 29
75 2 101 6 3 19 20 31
76 2 101 6 3 19 31 30
77 2 101 6 3 20 21 31
78 2 101 6 3 21 32 31
79 2 101 6 3 21 22 33
80 2 101 6 3 21 33 32
81 2 101 6 3 23 24 35
82 2 101 6 3 23 35 34
83 2 101 6 3 24 25 35
84 2 101 6 3 25 36 35
85 2 101 6 3 25 26 37
86 2 101 6 3 25 37 36
87 2 101 6 3 26 27 37
88 2 101 6 3 27 38 37
89 2 101 6 3 27 28 39
90 2 101 6 3 27 39 38
91 2 101 6 3 28 29 39
92 2 101 6 3 29 40 39
93 2 101 6 3 29 30 41
94 2 101 6 3 29 41 40
95 2 101 6 3 30 31 41
96 2 101 6 3 31 42 41
97 2 101 6 3 31 32 43
98 2 101 6 3 31 43 42
99 2 101 6 3 32 33 43
100 2 101 6 3 33 44 43
101 2 101 6 3 34 35 45
102 2 101 6 3 35 46 45
103 2 101 6 3 35 36 47
104 2 101 6 3 35 47 46
105 2 101 6 3 36 37 47
106 2 101 6 3 37 48 47
107 2 101 6 3 37 38 49
108 2 101 6 3 37 49 48
109 2 101 6 3 38 39 49
110 2 101 6 3 39 50 49
111 2 101 6 3 39 40 51
112 2 101 6 3 39 51 50
113 2 101 6 3 40 41 51
114 2 101 6 3 41 52 51
115 2 101 6 3 41 42 53
116 2 101 6 3 41 53 52
117 2 101 6 3 42 43 53
118 2 101 6 3 43 54 53
119 2 101 6 3 43 44 55
120 2 101 6 3 43 55 54
121 2 101 6 3 45 46 57
122 2 101 6 3 45 57 56
123 2 101 6 3 46 47 57
124 2 101 6 3 47 58 57
125 2 101 6 3 47 48 59
126 2 101 6 3 47 59 58
127 2 101 6 3 48 49 59
128 2 101 6 3 49 60 59
129 2 101 6 3 49 50 61
130 2 101 6 3 49 61 60
131 2 101 6 3 50 51 61
132 2 101 6 3 51 62 61
133 2 101 6 3 51 52 63
134 2 101 6 3 51 63 62
135 2 101 6 3 52 53 63
136 2 101 6 3 53 64 63
137 2 101 6 3 53 54 65
138 2 101 6 3 53 65 64
139 2 101 6 3 54 55 65
140 2 101 6 3 55 66 65
141 2 101 6 3 56 57 67
142 2 101 6 3 57 68 67
143 2 101 6 3 57 58 69
144 2 101 6 3 57 69 68
145 2 101 6 3 58 59 69
146 2 101 6 3 59 70 69
147 2 101 6 3 59 60 71
148 2 101 6 3 59 71 70
149 2 101 6 3 60 61 71
150 2 101 6 3 61 72 71
151 2 101 6 3 61 62 73
152 2 101 6 3 61 73 72
153 2 101 6 3 62 63 73
154 2 101 6 3 63 74 73
155 2 101 6 3 63 64 75
156 2 101 6 3 63 75 74
157 2 101 6 3 64 65 75
158 2 101 6 3 65 76 75
159 2 101 6 3 65 66 77
160 2 101 6 3 65 77 76
161 2 101 6 3 67 68 79
162 2 101 6 3 67 79 78
163 2 101 6 3 68 69 79
164 2 101 6 3 69 80 79
165 2 101 6 3 69 70 81
166 2 101 6 3 69 81 80
167 2 101 6 3 70 71 81
168 2 101 6 3 71 82 81
169 2 101 6 3 71 72 83
170 2 101 6 3 71 83 82
171 2 101 6 3 72 73 83
172 2 101 6 3 73 84 83
173 2 101 6 3 73 74 85
174 2 101 6 3 73 85 84
175 2 101 6 3 74 75 85
176 2 101 6 3 75 86 85
177 2 101 6 3 75 76 87
178 2 101 6 3 75 87 86
179 2 101 6 3 76 77 87
180 2 101 6 3 77 88 87
181 2 101 6 3 78 79 89
182 2 101 6 3 79 90 89
183 2 101 6 3 79 80 91
184 2 101 6 3 79 91 90
185 2 101 6 3 80 81 91
186 2 101 6 3 81 92 91
187 2 101 6 3 81 82 93
188 2 101 6 3 81 93 92
189 2 101 6 3 82 83 93
190 2 101 6 3 83 94 93
191 2 101 6 3 83 84 95
192 2 101 6 3 83 95 94
193 2 101 6 3 84 85 95
194 2 101 6 3 85 96 95
195 2 101 6 3 85 86 97
196 2 101 6 3 85 97 96
197 2 101 6 3 86 87 97
198 2 101 6 3 87 98 97
199 2 101 6 3 87 88 99
200 2 101 6 3 87 99 98
201 2 101 6 3 89 90 101
202 2 101 6 3 89 101 100
203 2 101 6 3 90 91 101
204 2 101 6 3 91 102 101
205 2 101 6 3 91 92 103
206 2 101 6 3 91 103 102
207 2 101 6 3 92 93 103
208 2 101 6 3 93 104 103
209 2 101 6 3 93 94 105
210 2 101 6 3 93 105 104
211 2 101 6 3 94 95 105
212 2 101 6 3 95 106 105
213 2 101 6 3 95 96 107
214 2 101 6 3 95 107 106
215 2 101 6 3 96 97 107
216 2 101 6 3 97 108 107
217 2 101 6 3 97 98 109
218 2 101 6 3 97 109 108
219 2 101 6 3 98 99 109
220 2 101 6 3 99 110 109
221 2 101 6 3 100 101 111
222 2 101 6 3 101 112 111
223 2 101 6 3 101 102 113
224 2 101 6 3 101 113 112
225 2 101 6 3 102 103 113
226 2 101 6 3 103 114 113
227 2 101 6 3 103 104 115
228 2 101 6 3 103 115 114
229 2 101 6 3 104 105 115
230 2 101 6 3 105 116 115
231 2 101 6 3 105 106 117
232 2 101 6 3 105 117 116
233 2 101 6 3 106 107 117
234 2 101 6 3 107 118 117
235 2 101 6 3 107 108 119
236 2 101 6 3 107 119 118
237 2 101 6 3 108 109 119
238 2 101 6 3 109 120 119
239 2 101 6 3 109 110 121
240 2 101 6 3 109 121 120
$ENDELM
#include "FastMarchingInterface.h"
#include "TLSAlgorithm.h"
#include "xAlgorithm.h"
#include "xLevelSet.h"
#include "xMesh.h"
#include "xPhysSurfByTagging.h"
#include "xRegion.h"
#include "xSimpleGeometry.h"
#include "xSpacePolynomial.h"
#include "xSubMesh.h"
#include "xDenseMatrix.h"
#include "xGenericSparseMatrix.h"
#include "xGenericSparseMatrixTraitPolicy.h"
namespace lalg {
typedef xGenericSparseMatrix<double,
xTraitMatrixNonSingular,
xTraitMatrixUnSym,
xTraitMatrixSparceCSC,
xTraitMatrixCindex> Matrix;
}
using namespace xfem;
using namespace lalg;
using namespace xtls;
using Trellis_Util::mPoint;
using AOMD::mEntity;
int main(int argc, char* argv[]) {
xMesh mesh(argv[1]);
// domain creation
xPlane p1(mPoint(0.1,0.1,0.), xVector(1.,0.,0.));
xPlane p2(mPoint(0.9,0.9,0.), xVector(-1.,0.,0.));
xLevelSet level_set(&mesh, xUnion(p1, p2));
xPhysSurfByTagging phys_surf(level_set);
xEntityFilter filter=std::bind1st(std::mem_fun(&xPhysSurfByTagging::strictOut), &phys_surf);
xSubMesh domain("domain", mesh);
xSubMesh bnd("bnd", mesh);
xSubMesh interior("interior", mesh);
xAcceptInSubMesh is_in_domain(domain);
xAcceptInSubMesh is_in_bnd(bnd);
for(xIter it=mesh.begin(2); it!=mesh.end(2); ++it) {
if(filter(*it)) {
domain.add(*it);
}
}
domain.modifyAllState();
for(xIter it=domain.begin(0); it!=domain.end(0); ++it) {
mEntity* n=*it;
for(int i=0; i<n->size(2); ++i) {
mEntity* e=n->get(2,i);
if(!is_in_domain(e)) {
bnd.add(n);
break;
}
}
}
for(xIter it=domain.begin(0); it!=domain.end(0); ++it) {
if(!is_in_bnd(*it)) {
interior.add(*it);
}
}
xRegion region(&domain);
// phi field
xDoubleManager double_manager;
xSpacePolynomialLagrange space("space", xSpace::SCALAR, 1);
xField field(&double_manager, space);
DeclareInterpolation(field, xValueCreator<xValueDouble>(), region.begin(), region.end());
DeclareState(field, xStateDofCreator<>(double_manager, "dofs"), region.begin(), region.end());
for(xIter it=region.begin(0); it!=region.end(0); ++it) {
xFiniteElementKeysOnly keys;
keys.setKeys(*it, field.begin(), field.end());
std::vector<xValue<double>*> vals;
vals.reserve(keys.sizeKey());
double_manager.getValPtr(keys.beginKey(), keys.endKey(), vals);
vals[0]->setVal(level_set(*it));
}
xEvalField<xIdentity<double> > eval_field(field);
xIntegrationRuleBasic integ_rule(1);
xExportGmshAscii pexport;
Export(eval_field, pexport, "phi_given", integ_rule, region.begin(), region.end());
std::function<double (const AOMD::mVertex&)> f_func=[](const AOMD::mVertex&){ return 1.; };
FastMarchingReinit(field,
region,
bnd.begin(0), bnd.end(0),
interior.begin(0), interior.end(0), f_func, false);
Export(eval_field, pexport, "phi", integ_rule, region.begin(), region.end());
// modes creation
xDenseMatrix coeffs;
FastMarchingModeExtension(field,
region,
bnd.begin(0), bnd.end(0),
interior.begin(0), interior.end(0),
coeffs);
// modes export
for(int i=0; i<bnd.size(0); ++i) {
int j=0;
for(xIter it=region.begin(0); it!=region.end(0); ++it, ++j) {
xFiniteElementKeysOnly keys;
keys.setKeys(*it, field.begin(), field.end());
std::vector<xValue<double>*> vals;
vals.reserve(keys.sizeKey());
double_manager.getValPtr(keys.beginKey(), keys.endKey(), vals);
vals[0]->setVal(coeffs(i,j));
}
std::ostringstream oss;
oss<<i;
Export(eval_field, pexport, "mode_"+oss.str(), integ_rule, region.begin(), region.end());
}
int j=0;
for(xIter it=region.begin(0); it!=region.end(0); ++it, ++j) {
xFiniteElementKeysOnly keys;
keys.setKeys(*it, field.begin(), field.end());
std::vector<xValue<double>*> vals;
vals.reserve(keys.sizeKey());
double_manager.getValPtr(keys.beginKey(), keys.endKey(), vals);
double sum=0.;
for(int i=0; i<bnd.size(0); ++i) {
sum+=coeffs(i,j);
}
vals[0]->setVal(sum);
}
Export(eval_field, pexport, "sum", integ_rule, region.begin(), region.end());
return 0;
}
/home/users/struct/fcazes/private/eXlibris/Applis/Plasticity/devel/patch_test/data
\ No newline at end of file
......@@ -12,10 +12,8 @@
#include "xDenseMatrix.h"
//
#ifdef HAVE_FASTMARCHING
// #include "gmshinputoutput.h"
#include "meshinterfacexRegion.h"
#include "FastMarching.h"
// #include "FM.h"
#endif
// std
#include <functional>
......@@ -83,23 +81,17 @@ void FastMarchingReinit(FIELD& phi,
}
meshinterfacexRegion mi(region);
entitystorage<meshinterfacexRegion, AOMD::mVertex, double> ls(mi);
// entitystorage<meshinterfacexRegion, AOMD::mVertex, vector3d<double> > gls(mi);
// entitystorage<meshinterfacexRegion, AOMD::mVertex, int> vn(mi);
for(ITERNODE it=known_begin; it!=known_end; ++it) {
auto v=*it;
std::vector<double> vals;
phi.getVals(v, vals);
ls.set(*static_cast<AOMD::mVertex*>(v), sign*vals[0]+max-min);
}
// for(ITERNODE it=tofind_begin; it!=tofind_end; ++it) {
// ls.set(*static_cast<AOMD::mVertex*>(*it), 1000.);
// }
entitytovertexiteratorconvertor<ITERNODE> known_conv_begin(known_begin);
entitytovertexiteratorconvertor<ITERNODE> known_conv_end(known_end);
entitytovertexiteratorconvertor<ITERNODE> tofind_conv_begin(tofind_begin);
entitytovertexiteratorconvertor<ITERNODE> tofind_conv_end(tofind_end);
fmeik(mi, ls, known_conv_begin, known_conv_end, tofind_conv_begin, tofind_conv_end, f_func);
// fmeik(mi, ls, known_conv_begin, known_conv_end, f_func, gls, vn);
for(ITERNODE it=tofind_begin; it!=tofind_end; ++it) {
auto v=*it;
double val;
......@@ -121,7 +113,6 @@ void FastMarchingModeExtension(FIELD& phi,
#else
meshinterfacexRegion mi(region);
entitystorage<meshinterfacexRegion, AOMD::mVertex, double> ls(mi);
// entitystorage<meshinterfacexRegion, AOMD::mVertex, vector3d<double> > gls(mi);
entitystorage<meshinterfacexRegion, AOMD::mVertex, modesValues> vn(mi);
const int nb_modes=std::distance(known_begin, known_end);
int i=0;
......@@ -134,18 +125,12 @@ void FastMarchingModeExtension(FIELD& phi,
mv[i]=1.;
vn.set(*static_cast<AOMD::mVertex*>(v), mv);
}
// for(ITERNODE it=tofind_begin; it!=tofind_end; ++it, ++i) {
// ls.set(*static_cast<AOMD::mVertex*>(*it), 1000.);
// modesValues mv(nb_modes, 0.);
// vn.set(*static_cast<AOMD::mVertex*>(*it), mv);
// }
entitytovertexiteratorconvertor<ITERNODE> known_conv_begin(known_begin);
entitytovertexiteratorconvertor<ITERNODE> known_conv_end(known_end);
entitytovertexiteratorconvertor<ITERNODE> tofind_conv_begin(tofind_begin);
entitytovertexiteratorconvertor<ITERNODE> tofind_conv_end(tofind_end);
std::function<double (const AOMD::mVertex&)> f_func=[](const AOMD::mVertex& v){return 1.;};
fmeik(mi, ls, known_conv_begin, known_conv_end, tofind_conv_begin, tofind_conv_end, f_func, vn);
// fmeik(mi, ls, known_conv_begin, known_conv_end, f_func, gls, vn);
coeffs.resize(nb_modes, region.size(0));
int j=0;
for(auto it=region.begin(0); it!=region.end(0); ++it, ++j) {
......
......@@ -2,7 +2,8 @@ MESH_FILE_TYPE = msh
PROCEDURE_PARAM_FILE = info.dat
MESH_FILE = mesh.msh
ZONE 10 = { MAT_CLASS = ElastoDam MAT_PARAM = mate.mat }
ZONE 101 = { MAT_CLASS = ElastoDam MAT_PARAM = mate.mat }
ZONE 102 = { MAT_CLASS = ElastoDam MAT_PARAM = mate.mat }
BC_LINE 6 = { DISPLACEMENT_Z FIX_AND_MEASURE = 1.e-3 }
BC_LINE 7 = { DISPLACEMENT_Z FIX = 0. }
BC_SURFACE 102 = { DISPLACEMENT_Z FIX_AND_MEASURE = 1.e-3 }
BC_LINE 12 = { DISPLACEMENT_Z FIX = 0. }
integ_order = 4
integ_order = 6
disp_space_type = PolyLagrange
disp_space_dim = V2Dxy
disp_space_order = 1
......@@ -10,7 +10,7 @@ delta_phi_max_ratio = 0.99
do_forced_delocalization = 0
do_crack_cut = 1
fit_to_vertices_ratio = -1.e-2
nb_step_max = 999
nb_step_max = 3000
export_manager = {
disp 10
eps_ref 10
......
......@@ -139,7 +139,11 @@ int main(int argc, char* argv[]) {
const double lambda=E*nu/(1.+nu)/(1.-2.*nu);
const double mu=E/2./(1.+nu);
std::cout<<"omega "<<omega<<" "<<omega*180/M_PI<<std::endl;
double alpha=0.5;
if(omega<208./180.*M_PI) {
alpha=0.75;
}
double f=sin(alpha*omega)+alpha*sin(omega);
for(int i=0; i<20 && fabs(f)>1.e-15; ++i) {
double df=omega*cos(alpha*omega)+sin(omega);
......@@ -150,7 +154,7 @@ int main(int argc, char* argv[]) {
post_pro.exportOnTime("singularity_exponent", 0., alpha);
const double coeff1=1./(((1.+alpha)*sin(omega*(1.+alpha)/2.))/((1.-alpha)*sin(omega*(1.-alpha)/2.))-1.)*1.e8;
const double coeff1=1./(((1.+alpha)*sin(omega*(1.+alpha)/2.))/((1.-alpha)*sin(omega*(1.-alpha)/2.))-1.)*1.e9;
const double coeff2=coeff1/(2.*mu);
EvalAsymptoticDisp eval_disp(omega, alpha, coeff2, lambda, mu);
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment