/[escript]/trunk/dudley/test/python/run_models.py
ViewVC logotype

Contents of /trunk/dudley/test/python/run_models.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 3774 - (show annotations)
Wed Jan 18 06:29:34 2012 UTC (7 years, 8 months ago) by jfenwick
File MIME type: text/x-python
File size: 77744 byte(s)
dudley, pasowrap, pycad

1 # -*- coding: utf-8 -*-
2
3 ########################################################
4 #
5 # Copyright (c) 2003-2010 by University of Queensland
6 # Earth Systems Science Computational Center (ESSCC)
7 # http://www.uq.edu.au/esscc
8 #
9 # Primary Business: Queensland, Australia
10 # Licensed under the Open Software License version 3.0
11 # http://www.opensource.org/licenses/osl-3.0.php
12 #
13 ########################################################
14
15 __copyright__="""Copyright (c) 2003-2010 by University of Queensland
16 Earth Systems Science Computational Center (ESSCC)
17 http://www.uq.edu.au/esscc
18 Primary Business: Queensland, Australia"""
19 __license__="""Licensed under the Open Software License version 3.0
20 http://www.opensource.org/licenses/osl-3.0.php"""
21 __url__="https://launchpad.net/escript-finley"
22
23 import unittest
24 import tempfile
25
26
27
28 VERBOSE=False and True
29
30 from esys.escript import *
31 from esys.escript.models import StokesProblemCartesian, PowerLaw, IncompressibleIsotropicFlowCartesian, FaultSystem, DarcyFlow
32 from esys.escript.models import Mountains
33 from esys.dudley import Rectangle, Brick
34
35 from math import pi
36 import numpy
37 import sys
38 import os
39 #====================================================================================================================
40 try:
41 DUDLEY_WORKDIR=os.environ['DUDLEY_WORKDIR']
42 except KeyError:
43 DUDLEY_WORKDIR='.'
44
45 #====================================================================================================================
46 class Test_StokesProblemCartesian2D(unittest.TestCase):
47 def setUp(self):
48 NE=6
49 self.TOL=1e-3
50 self.domain=Rectangle(NE,NE,order=-1)
51 def tearDown(self):
52 del self.domain
53 def test_PCG_P_0(self):
54 ETA=1.
55 P1=0.
56
57 x=self.domain.getX()
58 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
59 mask=whereZero(x[0]) * [1.,1.] \
60 +whereZero(x[0]-1) * [1.,1.] \
61 +whereZero(x[1]) * [1.,0.] \
62 +whereZero(x[1]-1) * [1.,1.]
63
64 sp=StokesProblemCartesian(self.domain)
65
66 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
67 u0=(1-x[0])*x[0]*[0.,1.]
68 p0=Scalar(-P1,ReducedSolution(self.domain))
69 sp.setTolerance(self.TOL)
70 u,p=sp.solve(u0*mask,p0,verbose=VERBOSE,max_iter=100,usePCG=True)
71
72 error_v0=Lsup(u[0]-u0[0])
73 error_v1=Lsup(u[1]-u0[1])/0.25
74 error_p=Lsup(p+P1*x[0]*x[1])
75 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
76 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
77 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
78
79 def test_PCG_P_small(self):
80 ETA=1.
81 P1=1.
82
83 x=self.domain.getX()
84 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
85 mask=whereZero(x[0]) * [1.,1.] \
86 +whereZero(x[0]-1) * [1.,1.] \
87 +whereZero(x[1]) * [1.,0.] \
88 +whereZero(x[1]-1) * [1.,1.]
89
90 sp=StokesProblemCartesian(self.domain)
91
92 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
93 u0=(1-x[0])*x[0]*[0.,1.]
94 p0=Scalar(-P1,ReducedSolution(self.domain))
95 sp.setTolerance(self.TOL)
96 u,p=sp.solve(u0,p0, verbose=VERBOSE,max_iter=100,usePCG=True)
97 error_v0=Lsup(u[0]-u0[0])
98 error_v1=Lsup(u[1]-u0[1])/0.25
99 error_p=Lsup(P1*x[0]*x[1]+p)
100 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
101 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
102 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
103
104 def test_PCG_P_large(self):
105 ETA=1.
106 P1=1000.
107
108 x=self.domain.getX()
109 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
110 mask=whereZero(x[0]) * [1.,1.] \
111 +whereZero(x[0]-1) * [1.,1.] \
112 +whereZero(x[1]) * [1.,0.] \
113 +whereZero(x[1]-1) * [1.,1.]
114
115 sp=StokesProblemCartesian(self.domain)
116
117 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
118 u0=(1-x[0])*x[0]*[0.,1.]
119 p0=Scalar(-P1,ReducedSolution(self.domain))
120 sp.setTolerance(self.TOL)
121 u,p=sp.solve(u0,p0, verbose=VERBOSE,max_iter=100,usePCG=True)
122 # u,p=sp.solve(u0,p0, verbose=VERBOSE,max_iter=100,usePCG=True)
123
124 error_v0=Lsup(u[0]-u0[0])
125 error_v1=Lsup(u[1]-u0[1])/0.25
126 error_p=Lsup(P1*x[0]*x[1]+p)/P1
127 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
128 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
129 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
130
131 def test_GMRES_P_0(self):
132 ETA=1.
133 P1=0.
134
135 x=self.domain.getX()
136 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
137 mask=whereZero(x[0]) * [1.,1.] \
138 +whereZero(x[0]-1) * [1.,1.] \
139 +whereZero(x[1]) * [1.,0.] \
140 +whereZero(x[1]-1) * [1.,1.]
141
142 sp=StokesProblemCartesian(self.domain)
143
144 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
145 u0=(1-x[0])*x[0]*[0.,1.]
146 p0=Scalar(-P1,ReducedSolution(self.domain))
147 sp.setTolerance(self.TOL)
148 u,p=sp.solve(u0,p0, verbose=VERBOSE,max_iter=50,usePCG=False,iter_restart=18)
149
150 error_v0=Lsup(u[0]-u0[0])
151 error_v1=Lsup(u[1]-u0[1])/0.25
152 error_p=Lsup(P1*x[0]*x[1]+p)
153 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
154 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
155 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
156
157 def test_GMRES_P_small(self):
158 ETA=1.
159 P1=1.
160
161 x=self.domain.getX()
162 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
163 mask=whereZero(x[0]) * [1.,1.] \
164 +whereZero(x[0]-1) * [1.,1.] \
165 +whereZero(x[1]) * [1.,0.] \
166 +whereZero(x[1]-1) * [1.,1.]
167
168 sp=StokesProblemCartesian(self.domain)
169
170 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
171 u0=(1-x[0])*x[0]*[0.,1.]
172 p0=Scalar(-P1,ReducedSolution(self.domain))
173 sp.setTolerance(self.TOL)
174 u,p=sp.solve(u0,p0, verbose=VERBOSE,max_iter=20,usePCG=False)
175
176 error_v0=Lsup(u[0]-u0[0])
177 error_v1=Lsup(u[1]-u0[1])/0.25
178 error_p=Lsup(P1*x[0]*x[1]+p)
179
180 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
181 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
182 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
183
184 def test_GMRES_P_large(self):
185 ETA=1.
186 P1=1000.
187
188 x=self.domain.getX()
189 F=-P1*x[1]*[1.,0]+(2*ETA-P1*x[0])*[0.,1.]
190 mask=whereZero(x[0]) * [1.,1.] \
191 +whereZero(x[0]-1) * [1.,1.] \
192 +whereZero(x[1]) * [1.,0.] \
193 +whereZero(x[1]-1) * [1.,1.]
194
195 sp=StokesProblemCartesian(self.domain)
196
197 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
198 u0=(1-x[0])*x[0]*[0.,1.]
199 p0=Scalar(-P1,ReducedSolution(self.domain))
200 sp.setTolerance(self.TOL)
201 u,p=sp.solve(u0,p0, verbose=VERBOSE,max_iter=100,usePCG=False)
202
203 error_v0=Lsup(u[0]-u0[0])
204 error_v1=Lsup(u[1]-u0[1])/0.25
205 error_p=Lsup(P1*x[0]*x[1]+p)/P1
206 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
207 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
208 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
209 #====================================================================================================================
210 class Test_StokesProblemCartesian3D(unittest.TestCase):
211 def setUp(self):
212 NE=6
213 self.TOL=1e-4
214 self.domain=Brick(NE,NE,NE,order=-1)
215 def tearDown(self):
216 del self.domain
217 def test_PCG_P_0(self):
218 ETA=1.
219 P1=0.
220
221 x=self.domain.getX()
222 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
223 x=self.domain.getX()
224 mask=whereZero(x[0]) * [1.,1.,1.] \
225 +whereZero(x[0]-1) * [1.,1.,1.] \
226 +whereZero(x[1]) * [1.,0.,1.] \
227 +whereZero(x[1]-1) * [1.,1.,1.] \
228 +whereZero(x[2]) * [1.,1.,0.] \
229 +whereZero(x[2]-1) * [1.,1.,1.]
230
231
232 sp=StokesProblemCartesian(self.domain)
233
234 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
235 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
236 p0=Scalar(-P1,ReducedSolution(self.domain))
237 sp.setTolerance(self.TOL)
238 u,p=sp.solve(u0,p0, verbose=VERBOSE ,max_iter=100,usePCG=True)
239
240 error_v0=Lsup(u[0]-u0[0])
241 error_v1=Lsup(u[1]-u0[1])
242 error_v2=Lsup(u[2]-u0[2])/0.25**2
243 error_p=Lsup(P1*x[0]*x[1]*x[2]+p)
244 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
245 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
246 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
247 self.assertTrue(error_v2<10*self.TOL, "2-velocity error too large.")
248
249 def test_PCG_P_small(self):
250 ETA=1.
251 P1=1.
252
253 x=self.domain.getX()
254 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
255 mask=whereZero(x[0]) * [1.,1.,1.] \
256 +whereZero(x[0]-1) * [1.,1.,1.] \
257 +whereZero(x[1]) * [1.,0.,1.] \
258 +whereZero(x[1]-1) * [1.,1.,1.] \
259 +whereZero(x[2]) * [1.,1.,0.] \
260 +whereZero(x[2]-1) * [1.,1.,1.]
261
262
263 sp=StokesProblemCartesian(self.domain)
264
265 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
266 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
267 p0=Scalar(-P1,ReducedSolution(self.domain))
268 sp.setTolerance(self.TOL)
269 u,p=sp.solve(u0,p0, verbose=VERBOSE ,max_iter=100,usePCG=True)
270 error_v0=Lsup(u[0]-u0[0])
271 error_v1=Lsup(u[1]-u0[1])
272 error_v2=Lsup(u[2]-u0[2])/0.25**2
273 error_p=Lsup(P1*x[0]*x[1]*x[2]+p)
274 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
275 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
276 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
277 self.assertTrue(error_v2<10*self.TOL, "2-velocity error too large.")
278
279 def test_PCG_P_large(self):
280 ETA=1.
281 P1=1000.
282
283 x=self.domain.getX()
284 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
285 mask=whereZero(x[0]) * [1.,1.,1.] \
286 +whereZero(x[0]-1) * [1.,1.,1.] \
287 +whereZero(x[1]) * [1.,0.,1.] \
288 +whereZero(x[1]-1) * [1.,1.,1.] \
289 +whereZero(x[2]) * [1.,1.,0.] \
290 +whereZero(x[2]-1) * [1.,1.,1.]
291
292
293 sp=StokesProblemCartesian(self.domain)
294
295 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
296 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
297 p0=Scalar(-P1,ReducedSolution(self.domain))
298 sp.setTolerance(self.TOL)
299 u,p=sp.solve(u0,p0, verbose=VERBOSE ,max_iter=100,usePCG=True)
300
301 error_v0=Lsup(u[0]-u0[0])
302 error_v1=Lsup(u[1]-u0[1])
303 error_v2=Lsup(u[2]-u0[2])/0.25**2
304 error_p=Lsup(P1*x[0]*x[1]*x[2]+p)/P1
305 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
306 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
307 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
308 self.assertTrue(error_v2<10*self.TOL, "2-velocity error too large.")
309
310 def test_GMRES_P_0(self):
311 ETA=1.
312 P1=0.
313
314 x=self.domain.getX()
315 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
316 x=self.domain.getX()
317 mask=whereZero(x[0]) * [1.,1.,1.] \
318 +whereZero(x[0]-1) * [1.,1.,1.] \
319 +whereZero(x[1]) * [1.,1.,1.] \
320 +whereZero(x[1]-1) * [1.,1.,1.] \
321 +whereZero(x[2]) * [1.,1.,0.] \
322 +whereZero(x[2]-1) * [1.,1.,1.]
323
324
325 sp=StokesProblemCartesian(self.domain)
326
327 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
328 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
329 p0=Scalar(-P1,ReducedSolution(self.domain))
330 sp.setTolerance(self.TOL)
331 u,p=sp.solve(u0,p0, verbose=VERBOSE,max_iter=100,usePCG=False,iter_restart=20)
332
333 error_v0=Lsup(u[0]-u0[0])
334 error_v1=Lsup(u[1]-u0[1])
335 error_v2=Lsup(u[2]-u0[2])/0.25**2
336 error_p=Lsup(P1*x[0]*x[1]*x[2]+p)
337 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
338 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
339 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
340 self.assertTrue(error_v2<10*self.TOL, "2-velocity error too large.")
341 def test_GMRES_P_small(self):
342 ETA=1.
343 P1=1.
344
345 x=self.domain.getX()
346 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
347 mask=whereZero(x[0]) * [1.,1.,1.] \
348 +whereZero(x[0]-1) * [1.,1.,1.] \
349 +whereZero(x[1]) * [1.,1.,1.] \
350 +whereZero(x[1]-1) * [1.,1.,1.] \
351 +whereZero(x[2]) * [1.,1.,0.] \
352 +whereZero(x[2]-1) * [1.,1.,1.]
353
354
355 sp=StokesProblemCartesian(self.domain)
356
357 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
358 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
359 p0=Scalar(-P1,ReducedSolution(self.domain))
360 sp.setTolerance(self.TOL/10)
361 u,p=sp.solve(u0,p0, verbose=VERBOSE,max_iter=100,usePCG=False)
362
363 error_v0=Lsup(u[0]-u0[0])
364 error_v1=Lsup(u[1]-u0[1])
365 error_v2=Lsup(u[2]-u0[2])/0.25**2
366 error_p=Lsup(P1*x[0]*x[1]*x[2]+p)/P1
367 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
368 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
369 self.assertTrue(error_v2<10*self.TOL, "2-velocity error too large.")
370 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
371 def test_GMRES_P_large(self):
372 ETA=1.
373 P1=1000.
374
375 x=self.domain.getX()
376 F=-P1*x[1]*x[2]*[1.,0.,0.]-P1*x[0]*x[2]*[0.,1.,0.]+(2*ETA*((1-x[0])*x[0]+(1-x[1])*x[1])-P1*x[0]*x[1])*[0.,0.,1.]
377 mask=whereZero(x[0]) * [1.,1.,1.] \
378 +whereZero(x[0]-1) * [1.,1.,1.] \
379 +whereZero(x[1]) * [1.,0.,1.] \
380 +whereZero(x[1]-1) * [1.,1.,1.] \
381 +whereZero(x[2]) * [1.,1.,0.] \
382 +whereZero(x[2]-1) * [1.,1.,1.]
383
384
385 sp=StokesProblemCartesian(self.domain)
386
387 sp.initialize(f=F,fixed_u_mask=mask,eta=ETA)
388 u0=(1-x[0])*x[0]*(1-x[1])*x[1]*[0.,0.,1.]
389 p0=Scalar(-P1,ReducedSolution(self.domain))
390 sp.setTolerance(self.TOL)
391 u,p=sp.solve(u0,p0, verbose=VERBOSE ,max_iter=100,usePCG=False)
392
393 error_v0=Lsup(u[0]-u0[0])
394 error_v1=Lsup(u[1]-u0[1])
395 error_v2=Lsup(u[2]-u0[2])/0.25**2
396 error_p=Lsup(P1*x[0]*x[1]*x[2]+p)/P1
397 self.assertTrue(error_p<10*self.TOL, "pressure error too large.")
398 self.assertTrue(error_v0<10*self.TOL, "0-velocity error too large.")
399 self.assertTrue(error_v1<10*self.TOL, "1-velocity error too large.")
400 self.assertTrue(error_v2<10*self.TOL, "2-velocity error too large.")
401 #====================================================================================================================
402 class Test_Darcy(unittest.TestCase):
403 # this is a simple test for the darcy flux problem
404 #
405 #
406 # p = 1/k * ( 1/2* (fb-f0)/xb* x **2 + f0 * x - ub*x ) + p0
407 #
408 # with f = (fb-f0)/xb* x + f0
409 #
410 # u = f - k * p,x = ub
411 #
412 # we prescribe pressure at x=x0=0 to p0
413 #
414 # if we prescribe the pressure on the bottom x=xb we set
415 #
416 # pb= 1/k * ( 1/2* (fb-f0)* xb + f0 * xb - ub*xb ) + p0 = 1/k * ((fb+f0)/2 - ub ) * xb + p0
417 #
418 # which leads to ub = (fb+f0)/2-k*(pb-p0)/xb
419 #
420 def rescaleDomain(self):
421 x=self.dom.getX().copy()
422 for i in range(self.dom.getDim()):
423 x_inf=inf(x[i])
424 x_sup=sup(x[i])
425 if i == self.dom.getDim()-1:
426 x[i]=-self.WIDTH*(x[i]-x_sup)/(x_inf-x_sup)
427 else:
428 x[i]=self.WIDTH*(x[i]-x_inf)/(x_sup-x_inf)
429 self.dom.setX(x)
430 def getScalarMask(self,include_bottom=True):
431 x=self.dom.getX().copy()
432 x_inf=inf(x[self.dom.getDim()-1])
433 x_sup=sup(x[self.dom.getDim()-1])
434 out=whereZero(x[self.dom.getDim()-1]-x_sup)
435 if include_bottom: out+=whereZero(x[self.dom.getDim()-1]-x_inf)
436 return wherePositive(out)
437 def getVectorMask(self,include_bottom=True):
438 x=self.dom.getX().copy()
439 out=Vector(0.,Solution(self.dom))
440 for i in range(self.dom.getDim()):
441 x_inf=inf(x[i])
442 x_sup=sup(x[i])
443 if i != self.dom.getDim()-1: out[i]+=whereZero(x[i]-x_sup)
444 if i != self.dom.getDim()-1 or include_bottom: out[i]+=whereZero(x[i]-x_inf)
445 return wherePositive(out)
446
447 def setSolutionFixedBottom(self, p0, pb, f0, fb, k):
448 d=self.dom.getDim()
449 x=self.dom.getX()[d-1]
450 xb=inf(x)
451 u=Vector(0.,Solution(self.dom))+kronecker(d)[d-1]*((f0+fb)/2.-k*(pb-p0)/xb)
452 p=1./k*((fb-f0)/(xb*2.)* x**2 - (fb-f0)/2.*x)+(pb-p0)/xb*x + p0
453 f= ((fb-f0)/xb* x + f0)*kronecker(Function(self.dom))[d-1]
454 return u,p,f
455
456 def testConstF_FixedBottom_smallK(self):
457 k=1.e-10
458 mp=self.getScalarMask(include_bottom=True)
459 mv=self.getVectorMask(include_bottom=False)
460 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=10.,k=k)
461 p=p_ref*mp
462 u=u_ref*mv
463 df=DarcyFlow(self.dom)
464 df.setValue(g=f,
465 location_of_fixed_pressure=mp,
466 location_of_fixed_flux=mv,
467 permeability=Scalar(k,Function(self.dom)))
468 df.setTolerance(rtol=self.TOL)
469 v,p=df.solve(u_ref,p, max_iter=100, verbose=VERBOSE)
470 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
471 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
472 def testConstF_FixedBottom_mediumK(self):
473 k=1.
474 mp=self.getScalarMask(include_bottom=True)
475 mv=self.getVectorMask(include_bottom=False)
476 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=10.,k=k)
477 p=p_ref*mp
478 u=u_ref*mv
479 df=DarcyFlow(self.dom)
480 df.setValue(g=f,
481 location_of_fixed_pressure=mp,
482 location_of_fixed_flux=mv,
483 permeability=Scalar(k,Function(self.dom)))
484 df.setTolerance(rtol=self.TOL)
485 v,p=df.solve(u,p,max_iter=100, verbose=VERBOSE )
486 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
487 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
488
489 def testConstF_FixedBottom_largeK(self):
490 k=1.e10
491 mp=self.getScalarMask(include_bottom=True)
492 mv=self.getVectorMask(include_bottom=False)
493 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=10.,k=k)
494 p=p_ref*mp
495 u=u_ref*mv
496 df=DarcyFlow(self.dom)
497 df.setValue(g=f,
498 location_of_fixed_pressure=mp,
499 location_of_fixed_flux=mv,
500 permeability=Scalar(k,Function(self.dom)))
501 df.setTolerance(rtol=self.TOL)
502 v,p=df.solve(u,p, max_iter=100, verbose=VERBOSE)
503 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
504 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
505
506 def testVarioF_FixedBottom_smallK(self):
507 k=1.e-10
508 mp=self.getScalarMask(include_bottom=True)
509 mv=self.getVectorMask(include_bottom=False)
510 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=30.,k=k)
511 p=p_ref*mp
512 u=u_ref*mv
513 df=DarcyFlow(self.dom)
514 #df.getSolverOptionsPressure().setVerbosityOn()
515 df.setValue(g=f,
516 location_of_fixed_pressure=mp,
517 location_of_fixed_flux=mv,
518 permeability=Scalar(k,Function(self.dom)))
519 df.setTolerance(rtol=self.TOL)
520 v,p=df.solve(u,p, max_iter=100, verbose=VERBOSE)
521
522 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
523 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
524
525 def testVarioF_FixedBottom_mediumK(self):
526 k=1.
527 mp=self.getScalarMask(include_bottom=True)
528 mv=self.getVectorMask(include_bottom=False)
529 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=30.,k=k)
530 p=p_ref*mp
531 u=u_ref*mv
532 df=DarcyFlow(self.dom)
533 df.setValue(g=f,
534 location_of_fixed_pressure=mp,
535 location_of_fixed_flux=mv,
536 permeability=Scalar(k,Function(self.dom)))
537 df.setTolerance(rtol=self.TOL)
538 v,p=df.solve(u,p, max_iter=100, verbose=VERBOSE)
539 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
540 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
541
542 def testVarioF_FixedBottom_largeK(self):
543 k=1.e10
544 mp=self.getScalarMask(include_bottom=True)
545 mv=self.getVectorMask(include_bottom=False)
546 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=30.,k=k)
547 p=p_ref*mp
548 u=u_ref*mv
549 df=DarcyFlow(self.dom)
550 df.setValue(g=f,
551 location_of_fixed_pressure=mp,
552 location_of_fixed_flux=mv,
553 permeability=Scalar(k,Function(self.dom)))
554 df.setTolerance(rtol=self.TOL)
555 v,p=df.solve(u,p, max_iter=100, verbose=VERBOSE)
556 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
557 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
558
559 def testConstF_FreeBottom_smallK(self):
560 k=1.e-10
561 mp=self.getScalarMask(include_bottom=False)
562 mv=self.getVectorMask(include_bottom=True)
563 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=10.,k=k)
564 p=p_ref*mp
565 u=u_ref*mv
566 df=DarcyFlow(self.dom)
567 df.setValue(g=f,
568 location_of_fixed_pressure=mp,
569 location_of_fixed_flux=mv,
570 permeability=Scalar(k,Function(self.dom)))
571 df.setTolerance(rtol=self.TOL)
572 v,p=df.solve(u,p, max_iter=100, verbose=VERBOSE)
573
574
575 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
576 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
577
578 def testConstF_FreeBottom_mediumK(self):
579 k=1.
580 mp=self.getScalarMask(include_bottom=False)
581 mv=self.getVectorMask(include_bottom=True)
582 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=10.,k=k)
583 p=p_ref*mp
584 u=u_ref*mv
585 df=DarcyFlow(self.dom)
586 df.setValue(g=f,
587 location_of_fixed_pressure=mp,
588 location_of_fixed_flux=mv,
589 permeability=Scalar(k,Function(self.dom)))
590 df.setTolerance(rtol=self.TOL)
591 v,p=df.solve(u,p, max_iter=100, verbose=VERBOSE)
592 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
593 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
594
595 def testConstF_FreeBottom_largeK(self):
596 k=1.e10
597 mp=self.getScalarMask(include_bottom=False)
598 mv=self.getVectorMask(include_bottom=True)
599 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=10.,k=k)
600 p=p_ref*mp
601 u=u_ref*mv
602 df=DarcyFlow(self.dom)
603 df.setValue(g=f,
604 location_of_fixed_pressure=mp,
605 location_of_fixed_flux=mv,
606 permeability=Scalar(k,Function(self.dom)))
607 df.setTolerance(rtol=self.TOL)
608 v,p=df.solve(u,p, max_iter=100, verbose=VERBOSE)
609 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
610 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
611
612 def testVarioF_FreeBottom_smallK(self):
613 k=1.e-10
614 mp=self.getScalarMask(include_bottom=False)
615 mv=self.getVectorMask(include_bottom=True)
616 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=30.,k=k)
617 p=p_ref*mp
618 u=u_ref*mv
619 df=DarcyFlow(self.dom)
620 df.setValue(g=f,
621 location_of_fixed_pressure=mp,
622 location_of_fixed_flux=mv,
623 permeability=Scalar(k,Function(self.dom)))
624 df.setTolerance(rtol=self.TOL)
625 v,p=df.solve(u,p, max_iter=100, verbose=VERBOSE)
626 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
627 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
628
629 def testVarioF_FreeBottom_mediumK(self):
630 k=1.
631 mp=self.getScalarMask(include_bottom=False)
632 mv=self.getVectorMask(include_bottom=True)
633 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=30.,k=k)
634 p=p_ref*mp
635 u=u_ref*mv
636 df=DarcyFlow(self.dom)
637 df.setValue(g=f,
638 location_of_fixed_pressure=mp,
639 location_of_fixed_flux=mv,
640 permeability=Scalar(k,Function(self.dom)))
641 df.setTolerance(rtol=self.TOL)
642 v,p=df.solve(u,p, max_iter=100, verbose=VERBOSE)
643 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
644 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
645
646 def testVarioF_FreeBottom_largeK(self):
647 k=1.e10
648 mp=self.getScalarMask(include_bottom=False)
649 mv=self.getVectorMask(include_bottom=True)
650 u_ref,p_ref,f=self.setSolutionFixedBottom(p0=2500,pb=4000.,f0=10.,fb=30.,k=k)
651 p=p_ref*mp
652 u=u_ref*mv
653 df=DarcyFlow(self.dom)
654 df.setValue(g=f,
655 location_of_fixed_pressure=mp,
656 location_of_fixed_flux=mv,
657 permeability=Scalar(k,Function(self.dom)))
658 df.setTolerance(rtol=self.TOL)
659 v,p=df.solve(u,p, max_iter=100, verbose=VERBOSE)
660 self.assertTrue(Lsup(v-u_ref)<self.TEST_TOL*Lsup(u_ref), "flux error too big.")
661 self.assertTrue(Lsup(p-p_ref)<self.TEST_TOL*Lsup(p_ref), "pressure error too big.")
662
663 class Test_Darcy2D(Test_Darcy):
664 TOL=1e-6
665 TEST_TOL=2.e-3
666 WIDTH=1.
667 def setUp(self):
668 NE=40 # wrning smaller NE may case a failure for VarioF tests due to discretization errors.
669 self.dom = Rectangle(NE/2,NE)
670 self.rescaleDomain()
671 def tearDown(self):
672 del self.dom
673 class Test_Darcy3D(Test_Darcy):
674 TOL=1e-6
675 WIDTH=1.
676 TEST_TOL=4.e-3
677 def setUp(self):
678 NE=25 # wrning smaller NE may case a failure for VarioF tests due to discretization errors.
679 self.dom = Brick(NE,NE,NE)
680 self.rescaleDomain()
681 def tearDown(self):
682 del self.dom
683
684 class Test_Rheologies(unittest.TestCase):
685 """
686 this is the program used to generate the powerlaw tests:
687
688 TAU_Y=100.
689 N=10
690 M=5
691
692 def getE(tau):
693 if tau<=TAU_Y:
694 return 1./(0.5+20*sqrt(tau))
695 else:
696 raise ValueError,"out of range."
697 tau=[ i* (TAU_Y/N) for i in range(N+1)] + [TAU_Y for i in range(M)]
698 e=[ tau[i]/getE(tau[i]) for i in range(N+1)]
699 e+= [ (i+1+j)* (max(e)/N) for j in range(M) ]
700
701 print tau
702 print e
703 """
704 TOL=1.e-8
705 def test_PowerLaw_Init(self):
706 pl=PowerLaw(numMaterials=3,verbose=VERBOSE)
707
708 self.assertEqual(pl.getNumMaterials(),3,"num materials is wrong")
709 self.assertEqual(pl.validMaterialId(0),True,"material id 0 not found")
710 self.assertEqual(pl.validMaterialId(1),True,"material id 1 not found")
711 self.assertEqual(pl.validMaterialId(2),True,"material id 2 not found")
712 self.assertEqual(pl.validMaterialId(3),False,"material id 3 not found")
713
714 self.assertEqual(pl.getFriction(),None,"initial friction wrong.")
715 self.assertEqual(pl.getTauY(),None,"initial tau y wrong.")
716 pl.setDruckerPragerLaw(tau_Y=10,friction=3)
717 self.assertEqual(pl.getFriction(),3,"friction wrong.")
718 self.assertEqual(pl.getTauY(),10,"tau y wrong.")
719
720 self.assertEqual(pl.getElasticShearModulus(),None,"initial shear modulus is wrong")
721 pl.setElasticShearModulus(1000)
722 self.assertEqual(pl.getElasticShearModulus(),1000.,"shear modulus is wrong")
723
724 e=pl.getEtaN()
725 self.assertEqual(len(e),3,"initial length of etaN is wrong.")
726 self.assertEqual(e,[None, None, None],"initial etaN are wrong.")
727 self.assertEqual(pl.getEtaN(0),None,"initial material id 0 is wrong")
728 self.assertEqual(pl.getEtaN(1),None,"initial material id 1 is wrong")
729 self.assertEqual(pl.getEtaN(2),None,"initial material id 2 is wrong")
730 self.assertRaises(ValueError, pl.getEtaN, 3)
731
732 self.assertRaises(ValueError,pl.setPowerLaws,eta_N=[10,20,30,40],tau_t=[2], power=[5,6,7])
733 self.assertRaises(ValueError,pl.setPowerLaws,eta_N=[10,20,30],tau_t=[2,4,5,5], power=[5,6,7])
734 self.assertRaises(ValueError,pl.setPowerLaws,eta_N=[10,20,30],tau_t=[2,1,2], power=[5,6,7,7])
735
736 pl.setPowerLaws(eta_N=[10,20,30],tau_t=[100,200,300], power=[1,2,3])
737 self.assertEqual(pl.getPower(),[1,2,3],"powers are wrong.")
738 self.assertEqual(pl.getTauT(),[100,200,300],"tau t are wrong.")
739 self.assertEqual(pl.getEtaN(),[10,20,30],"etaN are wrong.")
740
741 pl.setPowerLaw(id=0,eta_N=40,tau_t=400, power=4)
742 self.assertEqual(pl.getPower(),[4,2,3],"powers are wrong.")
743 self.assertEqual(pl.getTauT(),[400,200,300],"tau t are wrong.")
744 self.assertEqual(pl.getEtaN(),[40,20,30],"etaN are wrong.")
745
746 self.assertRaises(ValueError,pl.getPower,-1)
747 self.assertEqual(pl.getPower(0),4,"power 0 is wrong.")
748 self.assertEqual(pl.getPower(1),2,"power 1 is wrong.")
749 self.assertEqual(pl.getPower(2),3,"power 2 is wrong.")
750 self.assertRaises(ValueError,pl.getPower,3)
751
752 self.assertRaises(ValueError,pl.getTauT,-1)
753 self.assertEqual(pl.getTauT(0),400,"tau t 0 is wrong.")
754 self.assertEqual(pl.getTauT(1),200,"tau t 1 is wrong.")
755 self.assertEqual(pl.getTauT(2),300,"tau t 2 is wrong.")
756 self.assertRaises(ValueError,pl.getTauT,3)
757
758 self.assertRaises(ValueError,pl.getEtaN,-1)
759 self.assertEqual(pl.getEtaN(0),40,"eta n 0 is wrong.")
760 self.assertEqual(pl.getEtaN(1),20,"eta n 1 is wrong.")
761 self.assertEqual(pl.getEtaN(2),30,"eta n 2 is wrong.")
762 self.assertRaises(ValueError,pl.getEtaN,3)
763
764 def checkResult(self,id,gamma_dot_, eta, tau_ref):
765 self.assertTrue(eta>=0,"eta needs to be positive (test %s)"%id)
766 error=abs(gamma_dot_*eta-tau_ref)
767 self.assertTrue(error<=self.TOL*tau_ref,"eta is wrong: error = gamma_dot_*eta-tau_ref = %s * %s - %s = %s (test %s)"%(gamma_dot_,eta,tau_ref,error,id))
768
769 def test_PowerLaw_Linear(self):
770 taus= [0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 100.0, 100.0, 100.0, 100.0, 100.0]
771 gamma_dot_s=[0.0, 5.0, 10.0, 15.0, 20.0, 25.0, 30.0, 35.0, 40.0, 45.0, 50.0, 55.0, 60.0, 65.0, 70.0, 75.0]
772 pl=PowerLaw(numMaterials=1,verbose=VERBOSE)
773 pl.setDruckerPragerLaw(tau_Y=100.)
774 pl.setPowerLaw(eta_N=2.)
775 pl.setEtaTolerance(self.TOL)
776 for i in range(len(taus)): self.checkResult(i,gamma_dot_s[i], pl.getEtaEff(gamma_dot_s[i]),taus[i])
777
778 def test_PowerLaw_QuadLarge(self):
779 taus=[0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 100.0, 100.0, 100.0, 100.0, 100.0]
780 gamma_dot_s=[0.0, 405.0, 1610.0, 3615.0, 6420.0, 10025.0, 14430.0, 19635.0, 25640.0, 32445.0, 40050.0, 44055.0, 48060.0, 52065.0, 56070.0, 60075.0]
781 pl=PowerLaw(numMaterials=2,verbose=VERBOSE)
782 pl.setDruckerPragerLaw(tau_Y=100.)
783 pl.setPowerLaws(eta_N=[2.,0.01],tau_t=[1, 25.], power=[1,2])
784 pl.setEtaTolerance(self.TOL)
785 for i in range(len(taus)): self.checkResult(i,gamma_dot_s[i], pl.getEtaEff(gamma_dot_s[i]),taus[i])
786
787 def test_PowerLaw_QuadSmall(self):
788 taus=[0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 100.0, 100.0, 100.0, 100.0, 100.0]
789 gamma_dot_s=[0.0, 5.4, 11.6, 18.6, 26.4, 35.0, 44.4, 54.6, 65.6, 77.4, 90.0, 99.0, 108.0, 117.0, 126.0, 135.0]
790 pl=PowerLaw(numMaterials=2,verbose=VERBOSE)
791 pl.setDruckerPragerLaw(tau_Y=100.)
792 pl.setPowerLaws(eta_N=[2.,10.],tau_t=[1, 25.], power=[1,2])
793 pl.setEtaTolerance(self.TOL)
794 for i in range(len(taus)): self.checkResult(i,gamma_dot_s[i], pl.getEtaEff(gamma_dot_s[i]),taus[i])
795
796 def test_PowerLaw_CubeLarge(self):
797 taus=[0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 100.0, 100.0, 100.0, 100.0, 100.0]
798 gamma_dot_s=[0.0, 8.90625, 41.25, 120.46875, 270.0, 513.28125, 873.75, 1374.84375, 2040.0, 2892.65625, 3956.25, 4351.875, 4747.5, 5143.125, 5538.75, 5934.375]
799 pl=PowerLaw(numMaterials=2,verbose=VERBOSE)
800 pl.setDruckerPragerLaw(tau_Y=100.)
801 pl.setPowerLaws(eta_N=[2.,1./16.],tau_t=[1, 64.], power=[1,3])
802 pl.setEtaTolerance(self.TOL)
803 for i in range(len(taus)): self.checkResult(i,gamma_dot_s[i], pl.getEtaEff(gamma_dot_s[i]),taus[i])
804
805 def test_PowerLaw_CubeSmall(self):
806 taus=[0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 100.0, 100.0, 100.0, 100.0, 100.0]
807 gamma_dot_s=[0.0, 5.0390625, 10.3125, 16.0546875, 22.5, 29.8828125, 38.4375, 48.3984375, 60.0, 73.4765625, 89.0625, 97.96875, 106.875, 115.78125, 124.6875, 133.59375]
808 pl=PowerLaw(numMaterials=2,verbose=VERBOSE)
809 pl.setDruckerPragerLaw(tau_Y=100.)
810 pl.setPowerLaws(eta_N=[2.,25./4.],tau_t=[1, 64.], power=[1,3])
811 pl.setEtaTolerance(self.TOL)
812 for i in range(len(taus)): self.checkResult(i,gamma_dot_s[i], pl.getEtaEff(gamma_dot_s[i]),taus[i])
813
814 def test_PowerLaw_QuadLarge_CubeLarge(self):
815 taus=[0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 100.0, 100.0, 100.0, 100.0, 100.0]
816 gamma_dot_s=[0.0, 408.90625, 1641.25, 3720.46875, 6670.0, 10513.28125, 15273.75, 20974.84375, 27640.000000000004, 35292.65625, 43956.25, 48351.875, 52747.5, 57143.125, 61538.75, 65934.375]
817
818 pl=PowerLaw(numMaterials=3,verbose=VERBOSE)
819 pl.setDruckerPragerLaw(tau_Y=100.)
820 pl.setPowerLaws(eta_N=[2.,0.01,1./16.],tau_t=[1, 25.,64.], power=[1,2,3])
821 pl.setEtaTolerance(self.TOL)
822 for i in range(len(taus)): self.checkResult(i,gamma_dot_s[i], pl.getEtaEff(gamma_dot_s[i]),taus[i])
823
824 def test_PowerLaw_QuadLarge_CubeSmall(self):
825 taus=[0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 100.0, 100.0, 100.0, 100.0, 100.0]
826 gamma_dot_s=[0.0, 405.0390625, 1610.3125, 3616.0546875, 6422.5, 10029.8828125, 14438.4375, 19648.3984375, 25660.0, 32473.4765625, 40089.0625, 44097.96875, 48106.875, 52115.78125, 56124.6875, 60133.59375]
827
828 pl=PowerLaw(numMaterials=3,verbose=VERBOSE)
829 pl.setDruckerPragerLaw(tau_Y=100.)
830 pl.setPowerLaws(eta_N=[2.,0.01,25./4.],tau_t=[1, 25.,64.], power=[1,2,3])
831 pl.setEtaTolerance(self.TOL)
832 for i in range(len(taus)): self.checkResult(i,gamma_dot_s[i], pl.getEtaEff(gamma_dot_s[i]),taus[i])
833
834 def test_PowerLaw_QuadSmall_CubeLarge(self):
835 taus=[0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 100.0, 100.0, 100.0, 100.0, 100.0]
836 gamma_dot_s=[0.0, 9.30625, 42.85, 124.06875, 276.4, 523.28125, 888.15, 1394.44375, 2065.6, 2925.05625, 3996.25, 4395.875, 4795.5, 5195.125, 5594.75, 5994.375]
837 pl=PowerLaw(numMaterials=3,verbose=VERBOSE)
838 pl.setDruckerPragerLaw(tau_Y=100.)
839 pl.setPowerLaws(eta_N=[2.,10.,1./16.],tau_t=[1, 25.,64.], power=[1,2,3])
840 pl.setEtaTolerance(self.TOL)
841 for i in range(len(taus)): self.checkResult(i,gamma_dot_s[i], pl.getEtaEff(gamma_dot_s[i]),taus[i])
842
843 def test_PowerLaw_QuadSmall_CubeSmall(self):
844 taus=[0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 100.0, 100.0, 100.0, 100.0, 100.0]
845 gamma_dot_s=[0.0, 5.4390625, 11.9125, 19.6546875, 28.9, 39.8828125, 52.8375, 67.9984375, 85.6, 105.8765625, 129.0625, 141.96875, 154.875, 167.78125, 180.6875, 193.59375]
846 pl=PowerLaw(numMaterials=3,verbose=VERBOSE)
847 pl.setDruckerPragerLaw(tau_Y=100.)
848 pl.setPowerLaws(eta_N=[2.,10.,25./4.],tau_t=[1, 25.,64.], power=[1,2,3])
849 pl.setEtaTolerance(self.TOL)
850 for i in range(len(taus)): self.checkResult(i,gamma_dot_s[i], pl.getEtaEff(gamma_dot_s[i]),taus[i])
851
852 def test_PowerLaw_withShear(self):
853 taus=[0.0, 10.0, 20.0, 30.0, 40.0, 50.0, 60.0, 70.0, 80.0, 90.0, 100.0, 100.0, 100.0, 100.0, 100.0, 100.0]
854 gamma_dot_s=[0.0, 15.0, 30.0, 45.0, 60.0, 75.0, 90.0, 105.0, 120.0, 135.0, 150.0, 165.0, 180.0, 195.0, 210.0, 225.0]
855 pl=PowerLaw(numMaterials=1,verbose=VERBOSE)
856 pl.setDruckerPragerLaw(tau_Y=100.)
857 pl.setPowerLaw(eta_N=2.)
858 pl.setElasticShearModulus(3.)
859 dt=1./3.
860 pl.setEtaTolerance(self.TOL)
861 self.assertRaises(ValueError, pl.getEtaEff,gamma_dot_s[0])
862 for i in range(len(taus)): self.checkResult(i,gamma_dot_s[i], pl.getEtaEff(gamma_dot_s[i],dt=dt),taus[i])
863
864
865 class Test_FaultSystem(unittest.TestCase):
866 EPS=1.e-8
867 NE=10
868 def test_Fault_MaxValue(self):
869 dom=Rectangle(2*self.NE,2*self.NE)
870 x=dom.getX()
871 f=FaultSystem(dim=2)
872 f.addFault(V0=[0.5,0.], strikes=[3.*pi/4], ls=[0.70710678118654757], tag=1)
873 f.addFault(V0=[1.,0.5], strikes=[pi, pi/2], ls=[0.5,0.5], tag=2)
874
875 u=x[0]*(1.-x[0])*(1-x[1])
876 t, loc=f.getMaxValue(u)
877 p=f.getParametrization(x,t)[0]
878 m, l=loc(u), loc(p)
879 self.assertTrue( m == 0.25, "wrong max value")
880 self.assertTrue( t == 1, "wrong max tag")
881 self.assertTrue( l == 0., "wrong max location")
882
883 u=x[1]*(1.-x[1])*(1-x[0])*x[0]
884 t, loc=f.getMaxValue(u)
885 p=f.getParametrization(x,t)[0]
886 m, l=loc(u), loc(p)
887 self.assertTrue( m == 0.0625, "wrong max value")
888 self.assertTrue( t == 2, "wrong max tag")
889 self.assertTrue( l == 0.5, "wrong max location")
890
891 u=x[0]*(1.-x[0])*x[1]
892 t, loc=f.getMaxValue(u)
893 p=f.getParametrization(x,t)[0]
894 m, l=loc(u), loc(p)
895 self.assertTrue( m == 0.25, "wrong max value")
896 self.assertTrue( t == 2, "wrong max tag")
897 self.assertTrue( l == 1.0, "wrong max location")
898
899 u=x[1]*(1.-x[1])*x[0]
900 t, loc=f.getMaxValue(u)
901 p=f.getParametrization(x,t)[0]
902 m, l=loc(u), loc(p)
903 self.assertTrue( m == 0.25, "wrong max value")
904 self.assertTrue( t == 2, "wrong max tag")
905 self.assertTrue( l == 0., "wrong max location")
906
907 u=x[1]*(1.-x[1])*(1.-x[0])
908 t, loc=f.getMaxValue(u)
909 p=f.getParametrization(x,t)[0]
910 m, l=loc(u), loc(p)
911 self.assertTrue( m == 0.25, "wrong max value")
912 self.assertTrue( t == 1, "wrong max tag")
913 self.assertTrue( abs(l-0.70710678118654) <= self.EPS, "wrong max location")
914 def test_Fault_MinValue(self):
915 dom=Rectangle(2*self.NE,2*self.NE)
916 x=dom.getX()
917 f=FaultSystem(dim=2)
918 f.addFault(V0=[0.5,0.], strikes=[3.*pi/4], ls=[0.70710678118654757], tag=1)
919 f.addFault(V0=[1.,0.5], strikes=[pi, pi/2], ls=[0.5,0.5], tag=2)
920
921 u=-x[0]*(1.-x[0])*(1-x[1])
922 t, loc=f.getMinValue(u)
923 p=f.getParametrization(x,t)[0]
924 m, l=loc(u), loc(p)
925 self.assertTrue( m == -0.25, "wrong min value")
926 self.assertTrue( t == 1, "wrong min tag")
927 self.assertTrue( l == 0., "wrong min location")
928 u=-x[1]*(1.-x[1])*(1-x[0])*x[0]
929 t, loc=f.getMinValue(u)
930 p=f.getParametrization(x,t)[0]
931 m, l=loc(u), loc(p)
932 self.assertTrue( m == -0.0625, "wrong min value")
933 self.assertTrue( t == 2, "wrong min tag")
934 self.assertTrue( l == 0.5, "wrong min location")
935 u=-x[0]*(1.-x[0])*x[1]
936 t, loc=f.getMinValue(u)
937 p=f.getParametrization(x,t)[0]
938 m, l=loc(u), loc(p)
939 self.assertTrue( m == -0.25, "wrong min value")
940 self.assertTrue( t == 2, "wrong min tag")
941 self.assertTrue( l == 1.0, "wrong min location")
942 u=-x[1]*(1.-x[1])*x[0]
943 t, loc=f.getMinValue(u)
944 p=f.getParametrization(x,t)[0]
945 m, l=loc(u), loc(p)
946 self.assertTrue( m == -0.25, "wrong min value")
947 self.assertTrue( t == 2, "wrong min tag")
948 self.assertTrue( l == 0., "wrong min location")
949 u=-x[1]*(1.-x[1])*(1.-x[0])
950 t, loc=f.getMinValue(u)
951 p=f.getParametrization(x,t)[0]
952 m, l=loc(u), loc(p)
953 self.assertTrue( m == -0.25, "wrong min value")
954 self.assertTrue( t == 1, "wrong min tag")
955 self.assertTrue( abs(l-0.70710678118654) <= self.EPS, "wrong min location")
956
957
958 def test_Fault2D(self):
959 f=FaultSystem(dim=2)
960 top1=[ [1.,0.], [1.,1.], [0.,1.] ]
961 self.assertRaises(ValueError,f.addFault,V0=[1.,0],strikes=[pi/2, pi/2],ls=[1.,1.],tag=1,dips=top1)
962 f.addFault(V0=[1.,0],strikes=[pi/2, pi],ls=[1.,1.],tag=1)
963 self.assertTrue(f.getDim() == 2, "wrong dimension")
964 self.assertTrue( [ 1 ] == f.getTags(), "tags wrong")
965 self.assertTrue( 2. == f.getTotalLength(1), "length wrong")
966 self.assertTrue( 0. == f.getMediumDepth(1), "depth wrong")
967 self.assertTrue( (0., 2.) == f.getW0Range(1)," wrong W0 range")
968 self.assertTrue( (0., 0.) == f.getW1Range(1)," wrong W1 range")
969 self.assertTrue( [0., 1., 2.] == f.getW0Offsets(1)," wrong W0 offsets")
970 segs=f.getTopPolyline(1)
971 self.assertTrue( len(segs) == 3, "wrong number of segments")
972 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0")
973 self.assertTrue( segs[0].size == 2, "seg 0 has wrong size.")
974 self.assertTrue( numpy.linalg.norm(segs[0]-[1.,0.]) < self.EPS, "wrong vertex. 0 ")
975 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex 1")
976 self.assertTrue( segs[1].size == 2, "seg 1 has wrong size.")
977 self.assertTrue( numpy.linalg.norm(segs[1]-[1.,1.]) < self.EPS, "wrong vertex. 1 ")
978 self.assertTrue( isinstance(segs[2], numpy.ndarray), "wrong class of vertex 2")
979 self.assertTrue( segs[2].size == 2, "seg 2 has wrong size.")
980 self.assertTrue( numpy.linalg.norm(segs[2]-[0.,1.]) < self.EPS, "wrong vertex. 2 ")
981 c=f.getCenterOnSurface()
982 self.assertTrue( isinstance(c, numpy.ndarray), "center has wrong class")
983 self.assertTrue( c.size == 2, "center size is wrong")
984 self.assertTrue( numpy.linalg.norm(c-[2./3.,2./3.]) < self.EPS, "center has wrong coordinates.")
985 o=f.getOrientationOnSurface()/pi*180.
986 self.assertTrue( abs(o+45.) < self.EPS, "wrong orientation.")
987
988 top2=[ [10.,0.], [0.,10.] ]
989 f.addFault(V0=[10.,0],strikes=[3.*pi/4],ls=[14.142135623730951], tag=2, w0_offsets=[0,20], w1_max=20)
990 self.assertTrue( [ 1, 2 ] == f.getTags(), "tags wrong")
991 self.assertTrue( abs(f.getTotalLength(2)-14.1421356237) < self.EPS * 14.1421356237, "wrong length")
992 self.assertTrue( 0. == f.getMediumDepth(2), "depth wrong")
993 self.assertTrue( (0., 20.) == f.getW0Range(2)," wrong W0 range")
994 self.assertTrue( (0., 0.) == f.getW1Range(2)," wrong W1 range")
995 self.assertTrue( [0., 20.] == f.getW0Offsets(2)," wrong W0 offsets")
996 segs=f.getTopPolyline(2)
997 self.assertTrue( len(segs) == 2, "wrong number of segments")
998 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0")
999 self.assertTrue( numpy.linalg.norm(segs[0]-[10.,0.]) < self.EPS, "wrong vertex. 0 ")
1000 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex 1")
1001 self.assertTrue( numpy.linalg.norm(segs[1]-[0.,10.]) < self.EPS, "wrong vertex. 1 ")
1002 c=f.getCenterOnSurface()
1003 self.assertTrue( isinstance(c, numpy.ndarray), "center has wrong class")
1004 self.assertTrue( c.size == 2, "center size is wrong")
1005 self.assertTrue( numpy.linalg.norm(c-[12./5.,12./5.]) < self.EPS, "center has wrong coordinates.")
1006 o=f.getOrientationOnSurface()/pi*180.
1007 self.assertTrue( abs(o+45.) < self.EPS, "wrong orientation.")
1008
1009 s,d=f.getSideAndDistance([0.,0.], tag=1)
1010 self.assertTrue( s<0, "wrong side.")
1011 self.assertTrue( abs(d-1.)<self.EPS, "wrong distance.")
1012 s,d=f.getSideAndDistance([0.,2.], tag=1)
1013 self.assertTrue( s>0, "wrong side.")
1014 self.assertTrue( abs(d-1.)<self.EPS, "wrong distance.")
1015 s,d=f.getSideAndDistance([1.,2.], tag=1)
1016 self.assertTrue( s>0, "wrong side.")
1017 self.assertTrue( abs(d-1.)<self.EPS, "wrong distance.")
1018 s,d=f.getSideAndDistance([2.,1.], tag=1)
1019 self.assertTrue( s>0, "wrong side.")
1020 self.assertTrue( abs(d-1.)<self.EPS, "wrong distance.")
1021 s,d=f.getSideAndDistance([2.,0.], tag=1)
1022 self.assertTrue( s>0, "wrong side.")
1023 self.assertTrue( abs(d-1.)<self.EPS, "wrong distance.")
1024 s,d=f.getSideAndDistance([0.,-1.], tag=1)
1025 self.assertTrue( s<0, "wrong side.")
1026 self.assertTrue( abs(d-1.41421356237)<self.EPS, "wrong distance.")
1027 s,d=f.getSideAndDistance([-1.,0], tag=1)
1028 self.assertTrue( s<0, "wrong side.")
1029 self.assertTrue( abs(d-1.41421356237)<self.EPS, "wrong distance.")
1030
1031
1032 f.transform(rot=-pi/2., shift=[-1.,-1.])
1033 self.assertTrue( [ 1, 2 ] == f.getTags(), "tags after transformation wrong")
1034 self.assertTrue( 2. == f.getTotalLength(1), "length after transformation wrong")
1035 self.assertTrue( 0. == f.getMediumDepth(1), "depth after transformation wrong")
1036 self.assertTrue( (0., 2.) == f.getW0Range(1)," wrong W0 after transformation range")
1037 self.assertTrue( (0., 0.) == f.getW1Range(1)," wrong W1 rangeafter transformation ")
1038 self.assertTrue( [0., 1., 2.] == f.getW0Offsets(1)," wrong W0 offsetsafter transformation ")
1039 segs=f.getTopPolyline(1)
1040 self.assertTrue( len(segs) == 3, "wrong number of segmentsafter transformation ")
1041 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0 after transformation")
1042 self.assertTrue( segs[0].size == 2, "seg 0 has wrong size after transformation.")
1043 self.assertTrue( numpy.linalg.norm(segs[0]-[-1.,0.]) < self.EPS, "wrong vertex. 0 after transformation")
1044 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex after transformation1")
1045 self.assertTrue( segs[1].size == 2, "seg 1 has wrong size after transformation.")
1046 self.assertTrue( numpy.linalg.norm(segs[1]-[0.,0.]) < self.EPS, "wrong vertex. after transformation1 ")
1047 self.assertTrue( isinstance(segs[2], numpy.ndarray), "wrong class of vertex after transformation2")
1048 self.assertTrue( segs[2].size == 2, "seg 2 has wrong size after transformation.")
1049 self.assertTrue( numpy.linalg.norm(segs[2]-[0., 1.]) < self.EPS, "wrong vertex after transformation. 2 ")
1050 self.assertTrue( abs(f.getTotalLength(2)-14.1421356237) < self.EPS * 14.1421356237, "wrong length after transformation")
1051 self.assertTrue( 0. == f.getMediumDepth(2), "depth wrong after transformation")
1052 self.assertTrue( (0., 20.) == f.getW0Range(2)," wrong W0 range after transformation")
1053 self.assertTrue( (0., 0.) == f.getW1Range(2)," wrong W1 range after transformation")
1054 self.assertTrue( [0., 20.] == f.getW0Offsets(2)," wrong W0 offsets after transformation")
1055 segs=f.getTopPolyline(2)
1056 self.assertTrue( len(segs) == 2, "wrong number of segments after transformation")
1057 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0 after transformation")
1058 self.assertTrue( numpy.linalg.norm(segs[0]-[-1.,-9]) < self.EPS, "wrong vertex. 0 after transformation")
1059 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex 1 after transformation")
1060 self.assertTrue( numpy.linalg.norm(segs[1]-[9.,1.]) < self.EPS, "wrong vertex. 1 after transformation")
1061
1062 c=f.getCenterOnSurface()
1063 self.assertTrue( isinstance(c, numpy.ndarray), "center has wrong class")
1064 self.assertTrue( c.size == 2, "center size is wrong")
1065 self.assertTrue( numpy.linalg.norm(c-[7./5.,-7./5.]) < self.EPS, "center has wrong coordinates.")
1066 o=f.getOrientationOnSurface()/pi*180.
1067 self.assertTrue( abs(o-45.) < self.EPS, "wrong orientation.")
1068
1069 p=f.getParametrization([-1.,0.],1)
1070 self.assertTrue(p[1]==1., "wrong value.")
1071 self.assertTrue(abs(p[0])<self.EPS, "wrong value.")
1072 p=f.getParametrization([-0.5,0.],1)
1073 self.assertTrue(p[1]==1., "wrong value.")
1074 self.assertTrue(abs(p[0]-0.5)<self.EPS* 0.5, "wrong value.")
1075 p=f.getParametrization([0.,0.],1)
1076 self.assertTrue(p[1]==1., "wrong value.")
1077 self.assertTrue(abs(p[0]-1.)<self.EPS, "wrong value.")
1078 p=f.getParametrization([0.0000001,0.0000001],1, tol=1.e-8)
1079 self.assertTrue(p[1]==0., "wrong value.")
1080 p=f.getParametrization([0.0000001,0.0000001],1, tol=1.e-6)
1081 self.assertTrue(p[1]==1., "wrong value.")
1082 self.assertTrue(abs(p[0]-1.0000001)<self.EPS, "wrong value.")
1083 p=f.getParametrization([0.,0.5],1)
1084 self.assertTrue(p[1]==1., "wrong value.")
1085 self.assertTrue(abs(p[0]-1.5)<self.EPS, "wrong value.")
1086 p=f.getParametrization([0,1.],1)
1087 self.assertTrue(p[1]==1., "wrong value.")
1088 self.assertTrue(abs(p[0]-2.)<self.EPS, "wrong value.")
1089 p=f.getParametrization([1.,1.],1)
1090 self.assertTrue(p[1]==0., "wrong value.")
1091 p=f.getParametrization([0,1.11],1)
1092 self.assertTrue(p[1]==0., "wrong value.")
1093 p=f.getParametrization([-1,-9.],2)
1094 self.assertTrue(p[1]==1., "wrong value.")
1095 self.assertTrue(abs(p[0])<self.EPS, "wrong value.")
1096 p=f.getParametrization([9,1],2)
1097 self.assertTrue(p[1]==1., "wrong value.")
1098 self.assertTrue(abs(p[0]-20.)<self.EPS, "wrong value.")
1099
1100 def test_Fault3D(self):
1101 f=FaultSystem(dim=3)
1102 self.assertTrue(f.getDim() == 3, "wrong dimension")
1103
1104 top1=[ [0.,0.,0.], [1., 0., 0.] ]
1105 f.addFault(V0=[0.,0,0],strikes=[0.],ls=[1.,], dips=pi/2, depths=20.,tag=1)
1106 self.assertTrue( [ 1 ] == f.getTags(), "tags wrong")
1107 self.assertTrue( 1. == f.getTotalLength(1), "length wrong")
1108 self.assertTrue( 20. == f.getMediumDepth(1), "depth wrong")
1109 self.assertTrue( (0., 1.) == f.getW0Range(1)," wrong W0 range")
1110 self.assertTrue( (-20., 0.) == f.getW1Range(1)," wrong W1 range")
1111 self.assertTrue( [0., 1.] == f.getW0Offsets(1)," wrong W0 offsets")
1112 segs=f.getTopPolyline(1)
1113 self.assertTrue( len(segs) == 2, "wrong number of segments")
1114 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0")
1115 self.assertTrue( segs[0].size == 3, "seg 0 has wrong size.")
1116 self.assertTrue( numpy.linalg.norm(segs[0]-[0.,0.,0.]) < self.EPS, "wrong vertex. 0 ")
1117 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex 1")
1118 self.assertTrue( segs[1].size == 3, "seg 1 has wrong size.")
1119 self.assertTrue( numpy.linalg.norm(segs[1]-[1.,0.,0]) < self.EPS, "wrong vertex. 1 ")
1120 c=f.getCenterOnSurface()
1121 self.assertTrue( isinstance(c, numpy.ndarray), "center has wrong class")
1122 self.assertTrue( c.size == 3, "center size is wrong")
1123 self.assertTrue( numpy.linalg.norm(c-[0.5,0.,0.]) < self.EPS, "center has wrong coordinates.")
1124 o=f.getOrientationOnSurface()/pi*180.
1125 self.assertTrue( abs(o) < self.EPS, "wrong orientation.")
1126 d=f.getDips(1)
1127 self.assertTrue( len(d) == 1, "wrong number of dips")
1128 self.assertTrue( abs(d[0]-1.5707963267948966) < self.EPS, "wrong dip 0")
1129 sn=f.getSegmentNormals(1)
1130 self.assertTrue( len(sn) == 1, "wrong number of normals")
1131 self.assertTrue( isinstance(sn[0], numpy.ndarray), "wrong class of bottom vertex 0")
1132 self.assertTrue( numpy.linalg.norm(sn[0]-[0, -1., 0.]) < self.EPS, "wrong bottom vertex 1 ")
1133 dv=f.getDepthVectors(1)
1134 self.assertTrue( len(dv) == 2, "wrong number of depth vectors.")
1135 self.assertTrue( isinstance(dv[0], numpy.ndarray), "wrong class of depth vector 0")
1136 self.assertTrue( numpy.linalg.norm(dv[0]-[0., 0., -20.]) < self.EPS, "wrong depth vector 0 ")
1137 self.assertTrue( isinstance(dv[1], numpy.ndarray), "wrong class of depth vector 1")
1138 self.assertTrue( numpy.linalg.norm(dv[1]-[0., 0., -20.]) < self.EPS, "wrong depth vector 1 ")
1139 b=f.getBottomPolyline(1)
1140 self.assertTrue( len(b) == 2, "wrong number of bottom vertices")
1141 self.assertTrue( isinstance(b[0], numpy.ndarray), "wrong class of bottom vertex 0")
1142 self.assertTrue( numpy.linalg.norm(b[0]-[0., 0., -20.]) < self.EPS, "wrong bottom vertex 0 ")
1143 self.assertTrue( isinstance(b[1], numpy.ndarray), "wrong class of bottom vertex 1")
1144 self.assertTrue( numpy.linalg.norm(b[1]-[1., 0., -20.]) < self.EPS, "wrong bottom vertex 1 ")
1145 ds=f.getDepths(1)
1146 self.assertTrue( len(ds) == 2, "wrong number of depth")
1147 self.assertTrue( abs(ds[0]-20.) < self.EPS, "wrong depth at vertex 0 ")
1148 self.assertTrue( abs(ds[1]-20.) < self.EPS, "wrong depth at vertex 1 ")
1149
1150 top2=[ [0.,0.,0.], [0., 10., 0.] ]
1151 f.addFault(V0=[0.,0,0],strikes=[pi/2],ls=[10.,], dips=pi/2, depths=20.,tag=2)
1152 self.assertTrue( [ 1, 2 ] == f.getTags(), "tags wrong")
1153 self.assertTrue( 10. == f.getTotalLength(2), "length wrong")
1154 self.assertTrue( 20. == f.getMediumDepth(2), "depth wrong")
1155 self.assertTrue( (0., 10.) == f.getW0Range(2)," wrong W0 range")
1156 self.assertTrue( (-20., 0.) == f.getW1Range(2)," wrong W1 range")
1157 self.assertTrue( [0., 10.] == f.getW0Offsets(2)," wrong W0 offsets")
1158 segs=f.getTopPolyline(2)
1159 self.assertTrue( len(segs) == 2, "wrong number of segments")
1160 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0")
1161 self.assertTrue( segs[0].size == 3, "seg 0 has wrong size.")
1162 self.assertTrue( numpy.linalg.norm(segs[0]-[0.,0.,0.]) < self.EPS, "wrong vertex. 0 ")
1163 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex 1")
1164 self.assertTrue( segs[1].size == 3, "seg 1 has wrong size.")
1165 self.assertTrue( numpy.linalg.norm(segs[1]-[0., 10., 0]) < self.EPS, "wrong vertex. 1 ")
1166 d=f.getDips(2)
1167 self.assertTrue( len(d) == 1, "wrong number of dips")
1168 self.assertTrue( abs(d[0]-1.5707963267948966) < self.EPS, "wrong dip 0")
1169 sn=f.getSegmentNormals(2)
1170 self.assertTrue( len(sn) == 1, "wrong number of normals")
1171 self.assertTrue( isinstance(sn[0], numpy.ndarray), "wrong class of bottom vertex 0")
1172 self.assertTrue( numpy.linalg.norm(sn[0]-[1, 0., 0.]) < self.EPS, "wrong bottom vertex 1 ")
1173 dv=f.getDepthVectors(2)
1174 self.assertTrue( len(dv) == 2, "wrong number of depth vectors.")
1175 self.assertTrue( isinstance(dv[0], numpy.ndarray), "wrong class of depth vector 0")
1176 self.assertTrue( numpy.linalg.norm(dv[0]-[0., 0., -20.]) < self.EPS, "wrong depth vector 0 ")
1177 self.assertTrue( isinstance(dv[1], numpy.ndarray), "wrong class of depth vector 1")
1178 self.assertTrue( numpy.linalg.norm(dv[1]-[0., 0., -20.]) < self.EPS, "wrong depth vector 1 ")
1179 b=f.getBottomPolyline(2)
1180 self.assertTrue( len(b) == 2, "wrong number of bottom vertices")
1181 self.assertTrue( isinstance(b[0], numpy.ndarray), "wrong class of bottom vertex 0")
1182 self.assertTrue( numpy.linalg.norm(b[0]-[0., 0., -20.]) < self.EPS, "wrong bottom vertex 0 ")
1183 self.assertTrue( isinstance(b[1], numpy.ndarray), "wrong class of bottom vertex 1")
1184 self.assertTrue( numpy.linalg.norm(b[1]-[0., 10., -20.]) < self.EPS, "wrong bottom vertex 1 ")
1185 ds=f.getDepths(2)
1186 self.assertTrue( len(ds) == 2, "wrong number of depth")
1187 self.assertTrue( abs(ds[0]-20.) < self.EPS, "wrong depth at vertex 0 ")
1188 self.assertTrue( abs(ds[1]-20.) < self.EPS, "wrong depth at vertex 1 ")
1189
1190 top2=[ [10.,0.,0.], [0., 10., 0.] ]
1191 f.addFault(V0=[10.,0,0],strikes=3*pi/4,ls=14.142135623730951, dips=pi/2, depths=30.,tag=2)
1192 self.assertTrue( [ 1, 2 ] == f.getTags(), "tags wrong")
1193 self.assertTrue( abs(14.142135623730951 - f.getTotalLength(2)) <self.EPS, "length wrong")
1194 self.assertTrue( 30. == f.getMediumDepth(2), "depth wrong")
1195 self.assertTrue( (-30., 0.) == f.getW1Range(2)," wrong W1 range")
1196 segs=f.getTopPolyline(2)
1197 self.assertTrue( len(segs) == 2, "wrong number of segments")
1198 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0")
1199 self.assertTrue( segs[0].size == 3, "seg 0 has wrong size.")
1200 self.assertTrue( numpy.linalg.norm(segs[0]-[10.,0.,0.]) < self.EPS, "wrong vertex. 0 ")
1201 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex 1")
1202 self.assertTrue( segs[1].size == 3, "seg 1 has wrong size.")
1203 self.assertTrue( numpy.linalg.norm(segs[1]-[0., 10., 0]) < self.EPS, "wrong vertex. 1 ")
1204 d=f.getDips(2)
1205 self.assertTrue( len(d) == 1, "wrong number of dips")
1206 self.assertTrue( abs(d[0]-1.5707963267948966) < self.EPS, "wrong dip 0")
1207 sn=f.getSegmentNormals(2)
1208 self.assertTrue( len(sn) == 1, "wrong number of normals")
1209 self.assertTrue( isinstance(sn[0], numpy.ndarray), "wrong class of bottom vertex 0")
1210 self.assertTrue( numpy.linalg.norm(sn[0]-[0.70710678118654746, 0.70710678118654746, 0.]) < self.EPS, "wrong bottom vertex 1 ")
1211 dv=f.getDepthVectors(2)
1212 self.assertTrue( len(dv) == 2, "wrong number of depth vectors.")
1213 self.assertTrue( isinstance(dv[0], numpy.ndarray), "wrong class of depth vector 0")
1214 self.assertTrue( numpy.linalg.norm(dv[0]-[0., 0., -30.]) < self.EPS, "wrong depth vector 0 ")
1215 self.assertTrue( isinstance(dv[1], numpy.ndarray), "wrong class of depth vector 1")
1216 self.assertTrue( numpy.linalg.norm(dv[1]-[0., 0., -30.]) < self.EPS, "wrong depth vector 1 ")
1217 b=f.getBottomPolyline(2)
1218 self.assertTrue( len(b) == 2, "wrong number of bottom vertices")
1219 self.assertTrue( isinstance(b[0], numpy.ndarray), "wrong class of bottom vertex 0")
1220 self.assertTrue( numpy.linalg.norm(b[0]-[10., 0., -30.]) < self.EPS, "wrong bottom vertex 0 ")
1221 self.assertTrue( isinstance(b[1], numpy.ndarray), "wrong class of bottom vertex 1")
1222 self.assertTrue( numpy.linalg.norm(b[1]-[0., 10., -30.]) < self.EPS, "wrong bottom vertex 1 ")
1223 ds=f.getDepths(2)
1224 self.assertTrue( len(ds) == 2, "wrong number of depth")
1225 self.assertTrue( abs(ds[0]-30.) < self.EPS, "wrong depth at vertex 0 ")
1226 self.assertTrue( abs(ds[1]-30.) < self.EPS, "wrong depth at vertex 1 ")
1227
1228 top2=[ [10.,0.,0.], [0., 10., 0.] ]
1229 f.addFault(V0=[10.,0,0],strikes=3*pi/4,ls=14.142135623730951, dips=pi/4, depths=50.,tag=2)
1230 self.assertTrue( [ 1, 2 ] == f.getTags(), "tags wrong")
1231 self.assertTrue( abs(14.142135623730951 - f.getTotalLength(2)) <self.EPS, "length wrong")
1232 self.assertTrue( 50. == f.getMediumDepth(2), "depth wrong")
1233 self.assertTrue( (-50., 0.) == f.getW1Range(2)," wrong W1 range")
1234 segs=f.getTopPolyline(2)
1235 self.assertTrue( len(segs) == 2, "wrong number of segments")
1236 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0")
1237 self.assertTrue( segs[0].size == 3, "seg 0 has wrong size.")
1238 self.assertTrue( numpy.linalg.norm(segs[0]-[10.,0.,0.]) < self.EPS, "wrong vertex. 0 ")
1239 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex 1")
1240 self.assertTrue( segs[1].size == 3, "seg 1 has wrong size.")
1241 self.assertTrue( numpy.linalg.norm(segs[1]-[0., 10., 0]) < self.EPS, "wrong vertex. 1 ")
1242 d=f.getDips(2)
1243 self.assertTrue( len(d) == 1, "wrong number of dips")
1244 self.assertTrue( abs(d[0]-0.78539816339744828) < self.EPS, "wrong dip 0")
1245 sn=f.getSegmentNormals(2)
1246 self.assertTrue( len(sn) == 1, "wrong number of normals")
1247 self.assertTrue( isinstance(sn[0], numpy.ndarray), "wrong class of bottom vertex 0")
1248 self.assertTrue( numpy.linalg.norm(sn[0]-[0.5,0.5,0.70710678118654746]) < self.EPS, "wrong bottom vertex 1 ")
1249 dv=f.getDepthVectors(2)
1250 self.assertTrue( len(dv) == 2, "wrong number of depth vectors.")
1251 self.assertTrue( isinstance(dv[0], numpy.ndarray), "wrong class of depth vector 0")
1252 self.assertTrue( numpy.linalg.norm(dv[0]-[25., 25., -35.355339059327378]) < self.EPS, "wrong depth vector 0 ")
1253 self.assertTrue( isinstance(dv[1], numpy.ndarray), "wrong class of depth vector 1")
1254 self.assertTrue( numpy.linalg.norm(dv[1]-[25.,25., -35.355339059327378]) < self.EPS, "wrong depth vector 1 ")
1255 b=f.getBottomPolyline(2)
1256 self.assertTrue( len(b) == 2, "wrong number of bottom vertices")
1257 self.assertTrue( isinstance(b[0], numpy.ndarray), "wrong class of bottom vertex 0")
1258 self.assertTrue( numpy.linalg.norm(b[0]-[35., 25., -35.355339059327378]) < self.EPS, "wrong bottom vertex 0 ")
1259 self.assertTrue( isinstance(b[1], numpy.ndarray), "wrong class of bottom vertex 1")
1260 self.assertTrue( numpy.linalg.norm(b[1]-[25, 35., -35.355339059327378]) < self.EPS, "wrong bottom vertex 1 ")
1261 ds=f.getDepths(2)
1262 self.assertTrue( len(ds) == 2, "wrong number of depth")
1263 self.assertTrue( abs(ds[0]-50.) < self.EPS, "wrong depth at vertex 0 ")
1264 self.assertTrue( abs(ds[1]-50.) < self.EPS, "wrong depth at vertex 1 ")
1265
1266 top1=[ [10.,0.,0], [10.,10.,0], [0.,10.,0] ]
1267 f.addFault(V0=[10.,0.,0.],strikes=[pi/2, pi],ls=[10.,10.],tag=1, dips=pi/4, depths=20.)
1268 self.assertTrue( 20. == f.getTotalLength(1), "length wrong")
1269 self.assertTrue( 20. == f.getMediumDepth(1), "depth wrong")
1270 segs=f.getTopPolyline(1)
1271 self.assertTrue( len(segs) == 3, "wrong number of segments")
1272 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0")
1273 self.assertTrue( segs[0].size == 3, "seg 0 has wrong size.")
1274 self.assertTrue( numpy.linalg.norm(segs[0]-[10.,0.,0.]) < self.EPS, "wrong vertex. 0 ")
1275 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex 1")
1276 self.assertTrue( segs[1].size == 3, "seg 1 has wrong size.")
1277 self.assertTrue( numpy.linalg.norm(segs[1]-[10.,10.,0.]) < self.EPS, "wrong vertex. 1 ")
1278 self.assertTrue( isinstance(segs[2], numpy.ndarray), "wrong class of vertex 2")
1279 self.assertTrue( segs[2].size == 3, "seg 2 has wrong size.")
1280 self.assertTrue( numpy.linalg.norm(segs[2]-[0.,10.,0.]) < self.EPS, "wrong vertex. 2 ")
1281 d=f.getDips(1)
1282 self.assertTrue( len(d) == 2, "wrong number of dips")
1283 self.assertTrue( abs(d[0]-0.78539816339744828) < self.EPS, "wrong dip 0")
1284 self.assertTrue( abs(d[1]-0.78539816339744828) < self.EPS, "wrong dip 0")
1285 ds=f.getDepths(1)
1286 self.assertTrue( len(ds) == 3, "wrong number of depth")
1287 self.assertTrue( abs(ds[0]-20.) < self.EPS, "wrong depth at vertex 0 ")
1288 self.assertTrue( abs(ds[1]-20.) < self.EPS, "wrong depth at vertex 1 ")
1289 sn=f.getSegmentNormals(1)
1290 self.assertTrue( len(sn) == 2, "wrong number of normals")
1291 self.assertTrue( isinstance(sn[0], numpy.ndarray), "wrong class of bottom vertex 0")
1292 self.assertTrue( numpy.linalg.norm(sn[0]-[0.70710678118654746,0.,0.70710678118654746]) < self.EPS, "wrong bottom vertex 1 ")
1293 self.assertTrue( isinstance(sn[1], numpy.ndarray), "wrong class of bottom vertex 0")
1294 self.assertTrue( numpy.linalg.norm(sn[1]-[0.,0.70710678118654746,0.70710678118654746]) < self.EPS, "wrong bottom vertex 1 ")
1295 dv=f.getDepthVectors(1)
1296 self.assertTrue( len(dv) == 3, "wrong number of depth vectors.")
1297 self.assertTrue( isinstance(dv[0], numpy.ndarray), "wrong class of depth vector 0")
1298 self.assertTrue( numpy.linalg.norm(dv[0]-[14.142135623730951, 0., -14.142135623730951]) < self.EPS, "wrong depth vector 0 ")
1299 self.assertTrue( isinstance(dv[1], numpy.ndarray), "wrong class of depth vector 1")
1300 self.assertTrue( numpy.linalg.norm(dv[1]-[11.547005383792515,11.547005383792515, -11.547005383792515]) < self.EPS, "wrong depth vector 2 ")
1301 self.assertTrue( isinstance(dv[2], numpy.ndarray), "wrong class of depth vector 1")
1302 self.assertTrue( numpy.linalg.norm(dv[2]-[0.,14.142135623730951, -14.142135623730951]) < self.EPS, "wrong depth vector 2 ")
1303 segs=f.getBottomPolyline(1)
1304 self.assertTrue( len(segs) == 3, "wrong number of segments")
1305 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0")
1306 self.assertTrue( segs[0].size == 3, "seg 0 has wrong size.")
1307 self.assertTrue( numpy.linalg.norm(segs[0]-[24.142135623730951,0.,-14.142135623730951]) < self.EPS, "wrong vertex. 0 ")
1308 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex 1")
1309 self.assertTrue( segs[1].size == 3, "seg 1 has wrong size.")
1310 self.assertTrue( numpy.linalg.norm(segs[1]-[21.547005383792515,21.547005383792515, -11.547005383792515]) < self.EPS, "wrong vertex. 1 ")
1311 self.assertTrue( isinstance(segs[2], numpy.ndarray), "wrong class of vertex 2")
1312 self.assertTrue( segs[2].size == 3, "seg 2 has wrong size.")
1313 self.assertTrue( numpy.linalg.norm(segs[2]-[0., 24.142135623730951, -14.142135623730951]) < self.EPS, "wrong vertex. 2 ")
1314 self.assertTrue( abs(0.-f.getW0Range(1)[0]) <=self.EPS," wrong W0 range (0)")
1315 self.assertTrue( abs(31.857329272664341-f.getW0Range(1)[1]) <=self.EPS," wrong W0 range (1)")
1316 self.assertTrue( abs(-20. - f.getW1Range(1)[0]) <=self.EPS," wrong W1 range (0)")
1317 self.assertTrue( abs(0. - f.getW1Range(1)[1]) <=self.EPS," wrong W1 range (1)")
1318 self.assertTrue( abs(0.0-f.getW0Offsets(1)[0])<=self.EPS," wrong W0 offsets (0)")
1319 self.assertTrue( abs(15.92866463633217-f.getW0Offsets(1)[1])<=self.EPS," wrong W0 offsets (1)")
1320 self.assertTrue( abs(31.857329272664341-f.getW0Offsets(1)[2])<=self.EPS," wrong W0 offsets(2)")
1321 #
1322 # ============ fresh start ====================
1323 #
1324 f.addFault(V0=[1.,0,0.],strikes=[pi/2, pi],ls=[1.,1.], dips=pi/2,depths=20,tag=1)
1325 f.addFault(V0=[10.,0,0],strikes=[3.*pi/4],ls=[14.142135623730951], tag=2, w0_offsets=[0,20], w1_max=20, dips=pi/2,depths=20)
1326 c=f.getCenterOnSurface()
1327 self.assertTrue( isinstance(c, numpy.ndarray), "center has wrong class")
1328 self.assertTrue( c.size == 3, "center size is wrong")
1329 self.assertTrue( numpy.linalg.norm(c-[12./5.,12./5.,0.]) < self.EPS, "center has wrong coordinates.")
1330 o=f.getOrientationOnSurface()/pi*180.
1331 self.assertTrue( abs(o+45.) < self.EPS, "wrong orientation.")
1332
1333 f.transform(rot=-pi/2., shift=[-1.,-1.,0.])
1334 self.assertTrue( [ 1, 2 ] == f.getTags(), "tags after transformation wrong")
1335 self.assertTrue( 2. == f.getTotalLength(1), "length after transformation wrong")
1336 self.assertTrue( 20. == f.getMediumDepth(1), "depth after transformation wrong")
1337 rw0=f.getW0Range(1)
1338 self.assertTrue( len(rw0) ==2, "wo range has wrong length")
1339 self.assertTrue( abs(rw0[0]) < self.EPS,"W0 0 wrong.")
1340 self.assertTrue( abs(rw0[1]-2.) < self.EPS,"W0 1 wrong.")
1341 self.assertTrue( (-20., 0.) == f.getW1Range(1)," wrong W1 rangeafter transformation ")
1342 dips=f.getDips(1)
1343 self.assertTrue(len(dips) == 2, "wrong number of dips.")
1344 self.assertTrue( abs(dips[0]-1.5707963267948966) <= self.EPS, "wrong dip")
1345 self.assertTrue( abs(dips[1]-1.5707963267948966) <= self.EPS, "wrong dip")
1346 ds=f.getDepths(1)
1347 self.assertTrue( len(ds) == 3, "wrong number of depth")
1348 self.assertTrue( abs(ds[0]-20.) < self.EPS, "wrong depth at vertex 0 ")
1349 self.assertTrue( abs(ds[1]-20.) < self.EPS, "wrong depth at vertex 1 ")
1350 self.assertTrue( abs(ds[2]-20.) < self.EPS, "wrong depth at vertex 1 ")
1351 segs=f.getTopPolyline(1)
1352 self.assertTrue( len(segs) == 3, "wrong number of segmentsafter transformation ")
1353 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0 after transformation")
1354 self.assertTrue( segs[0].size == 3, "seg 0 has wrong size after transformation.")
1355 self.assertTrue( numpy.linalg.norm(segs[0]-[-1.,0.,0.]) < self.EPS, "wrong vertex. 0 after transformation")
1356 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex after transformation1")
1357 self.assertTrue( segs[1].size == 3, "seg 1 has wrong size after transformation.")
1358 self.assertTrue( numpy.linalg.norm(segs[1]-[0.,0.,0.]) < self.EPS, "wrong vertex. after transformation1 ")
1359 self.assertTrue( isinstance(segs[2], numpy.ndarray), "wrong class of vertex after transformation2")
1360 self.assertTrue( segs[2].size == 3, "seg 2 has wrong size after transformation.")
1361 self.assertTrue( numpy.linalg.norm(segs[2]-[0., 1.,0.]) < self.EPS, "wrong vertex after transformation. 2 ")
1362 self.assertTrue( abs(f.getTotalLength(2)-14.1421356237) < self.EPS * 14.1421356237, "wrong length after transformation")
1363 self.assertTrue( 20. == f.getMediumDepth(2), "depth wrong after transformation")
1364 rw0=f.getW0Range(2)
1365 self.assertTrue( len(rw0) ==2, "wo range has wrong length")
1366 self.assertTrue( abs(rw0[0]) < self.EPS,"W0 0 wrong.")
1367 self.assertTrue( abs(rw0[1]-20.) < self.EPS,"W0 1 wrong.")
1368 self.assertTrue( (-20., 0.) == f.getW1Range(2)," wrong W1 range after transformation")
1369 self.assertTrue( [0., 20.] == f.getW0Offsets(2)," wrong W0 offsets after transformation")
1370 dips=f.getDips(2)
1371 self.assertTrue(len(dips) == 1, "wrong number of dips.")
1372 self.assertTrue( abs(dips[0]-1.5707963267948966) <= self.EPS, "wrong dip")
1373 ds=f.getDepths(2)
1374 self.assertTrue( len(ds) == 2, "wrong number of depth")
1375 self.assertTrue( abs(ds[0]-20.) < self.EPS, "wrong depth at vertex 0 ")
1376 self.assertTrue( abs(ds[1]-20.) < self.EPS, "wrong depth at vertex 1 ")
1377 segs=f.getTopPolyline(2)
1378 self.assertTrue( len(segs) == 2, "wrong number of segments after transformation")
1379 self.assertTrue( isinstance(segs[0], numpy.ndarray), "wrong class of vertex 0 after transformation")
1380 self.assertTrue( numpy.linalg.norm(segs[0]-[-1.,-9,0.]) < self.EPS, "wrong vertex. 0 after transformation")
1381 self.assertTrue( isinstance(segs[1], numpy.ndarray), "wrong class of vertex 1 after transformation")
1382 self.assertTrue( numpy.linalg.norm(segs[1]-[9.,1.,0.]) < self.EPS, "wrong vertex. 1 after transformation")
1383 #
1384 # ============ fresh start ====================
1385 #
1386 f=FaultSystem(dim=3)
1387
1388 top1=[ [0.,0.,0.], [1., 0., 0.] ]
1389 f.addFault(V0=[0.,0,0],strikes=[0.],ls=[1.,], dips=pi/2, depths=1.,tag=1)
1390 top1=[ [10.,0.,0], [10.,10.,0], [0.,10.,0] ]
1391 f.addFault(V0=[10.,0.,0.],strikes=[pi/2, pi],ls=[10.,10.],tag=2, dips=pi/4, depths=20.)
1392
1393 p,m=f.getParametrization([0.3,0.,-0.5],1)
1394 self.assertTrue(length(p-[0.3,-0.5]) <= self.EPS, "wrong value.")
1395 self.assertTrue(m==1., "wrong value.")
1396
1397 p,m=f.getParametrization([0.5,0.,-0.5],1)
1398 self.assertTrue(length(p-[0.5,-0.5]) <= self.EPS, "wrong value.")
1399 self.assertTrue(m==1., "wrong value.")
1400
1401 p,m=f.getParametrization([0.25,0.,-0.5],1)
1402 self.assertTrue(length(p-[0.25,-0.5]) <= self.EPS, "wrong value.")
1403 self.assertTrue(m==1., "wrong value.")
1404
1405 p,m=f.getParametrization([0.5,0.,-0.25],1)
1406 self.assertTrue(length(p-[0.5,-0.25]) <= self.EPS, "wrong value.")
1407 self.assertTrue(m==1., "wrong value.")
1408
1409 p,m=f.getParametrization([0.001,0.,-0.001],1)
1410 self.assertTrue(length(p-[0.001, -0.001]) <= self.EPS, "wrong value.")
1411 self.assertTrue(m==1., "wrong value.")
1412
1413 p,m=f.getParametrization([0.001,0.,0.001],1)
1414 self.assertTrue(m==0., "wrong value.")
1415
1416 p,m=f.getParametrization([0.999,0.,0.001],1)
1417 self.assertTrue(m==0., "wrong value.")
1418
1419 p,m=f.getParametrization([1.001,0.,-0.001],1)
1420 self.assertTrue(m==0., "wrong value.")
1421 p,m=f.getParametrization([1.001,0.,-0.1],1)
1422 self.assertTrue(m==0., "wrong value.")
1423 p,m=f.getParametrization([1.001,0.,-0.000000001],1)
1424 self.assertTrue(m==0., "wrong value.")
1425
1426 p,m=f.getParametrization([0.999,0.,-0.001],1)
1427 self.assertTrue(length(p-[0.999, -0.001]) <= self.EPS, "wrong value.")
1428 self.assertTrue(m==1., "wrong value.")
1429
1430 p,m=f.getParametrization([ 16.29252873 , 6.46410161 ,-6.29252873], 2, tol=1.e-7)
1431 self.assertTrue(m==1., "wrong value.")
1432 self.assertTrue(length(p-[4.7785993908996511, -10]) <= self.EPS*10., "wrong value.")
1433 p,m=f.getParametrization([15.77350269, 12.77350269, -5.77350269], 2, tol=1.e-7)
1434 self.assertTrue(m==1., "wrong value.")
1435 self.assertTrue(length(p-[11.150065245432518, -10]) <= self.EPS*10., "wrong value.")
1436
1437 p,m=f.getParametrization([ 3., 17.0710678, -7.0710678], 2, tol=1.e-7)
1438 self.assertTrue(m==1., "wrong value.")
1439 self.assertTrue(length(p-[27.078729881764687, -10]) <= self.EPS*10., "wrong value.")
1440 p,m=f.getParametrization([9.30940108, 16.55204176, -6.55204176], 2, tol=1.e-7)
1441 self.assertTrue(m==1., "wrong value.")
1442 self.assertTrue(length(p-[20.707264027231822, -10]) <= self.EPS*10., "wrong value.")
1443
1444 p,m=f.getParametrization([ 21.54700538, 21.54700538, -11.54700538], 2, tol=1.e-7)
1445 self.assertTrue(m==1., "wrong value.")
1446 self.assertTrue(length(p-[15.92866463633217, -20]) <= self.EPS*10., "wrong value.")
1447
1448 p,m=f.getParametrization([ 0.,0.,0.], 2, tol=1.e-7)
1449 self.assertTrue(m==0., "wrong value.")
1450
1451 p,m=f.getParametrization([ 11.,11.,0.], 2, tol=1.e-7)
1452 self.assertTrue(m==0., "wrong value.")
1453
1454
1455 s,d=f.getSideAndDistance([0.,-1.,0.], tag=1)
1456 self.assertTrue( s>0, "wrong side.")
1457 self.assertTrue( abs(d-1.)<self.EPS, "wrong distance.")
1458 s,d=f.getSideAndDistance([1.,-1.,0.], tag=1)
1459 self.assertTrue( s>0, "wrong side.")
1460 self.assertTrue( abs(d-1.)<self.EPS, "wrong distance.")
1461 s,d=f.getSideAndDistance([0.,1.,0.], tag=1)
1462 self.assertTrue( s<0, "wrong side.")
1463 self.assertTrue( abs(d-1.)<self.EPS, "wrong distance.")
1464 s,d=f.getSideAndDistance([1.,1.,0.], tag=1)
1465 self.assertTrue( s<0, "wrong side.")
1466 self.assertTrue( abs(d-1.)<self.EPS, "wrong distance.")
1467
1468
1469 s,d=f.getSideAndDistance([0.,0.,0.], tag=2)
1470 self.assertTrue( s<0, "wrong side.")
1471 self.assertTrue( abs(d-10.)<self.EPS, "wrong distance.")
1472 s,d=f.getSideAndDistance([5.,5.,0.], tag=2)
1473 self.assertTrue( s<0, "wrong side.")
1474 self.assertTrue( abs(d-5.)<self.EPS, "wrong distance.")
1475
1476 s,d=f.getSideAndDistance([10.,10.,-1.], tag=2)
1477 self.assertTrue( s<0, "wrong side.")
1478 self.assertTrue( abs(d-0.70710678118654757)<self.EPS, "wrong distance.")
1479 s,d=f.getSideAndDistance([10.,10.,-2.], tag=2)
1480 self.assertTrue( s<0, "wrong side.")
1481 self.assertTrue( abs(d-2.*0.70710678118654757)<self.EPS, "wrong distance.")
1482 s,d=f.getSideAndDistance([10.,10.,-3.], tag=2)
1483 self.assertTrue( s<0, "wrong side.")
1484 self.assertTrue( abs(d-3.*0.70710678118654757)<self.EPS, "wrong distance.")
1485
1486 s,d=f.getSideAndDistance([5.,12.,0], tag=2)
1487 self.assertTrue( s>0, "wrong side.")
1488 self.assertTrue( abs(d-2*0.70710678118654757)<self.EPS, "wrong distance.")
1489 s,d=f.getSideAndDistance([5.,12.,-1], tag=2)
1490 self.assertTrue( s>0, "wrong side.")
1491 self.assertTrue( abs(d-0.70710678118654757)<self.EPS, "wrong distance.")
1492 s,d=f.getSideAndDistance([5.,12.,-2], tag=2)
1493 # s not checked as it is undefined.
1494 self.assertTrue( abs(d)<self.EPS, "wrong distance.")
1495 s,d=f.getSideAndDistance([5.,12.,-3], tag=2)
1496 self.assertTrue( s<0, "wrong side.")
1497 self.assertTrue( abs(d-0.70710678118654757)<self.EPS, "wrong distance.")
1498 s,d=f.getSideAndDistance([5.,12.,-4], tag=2)
1499 self.assertTrue( s<0, "wrong side.")
1500 self.assertTrue( abs(d-2.*0.70710678118654757)<self.EPS, "wrong distance.")
1501
1502 if __name__ == '__main__':
1503 suite = unittest.TestSuite()
1504 suite.addTest(unittest.makeSuite(Test_FaultSystem))
1505 # suite.addTest(unittest.makeSuite(Test_StokesProblemCartesian2D))
1506 suite.addTest(unittest.makeSuite(Test_Darcy3D))
1507 suite.addTest(unittest.makeSuite(Test_Darcy2D))
1508 # suite.addTest(Test_StokesProblemCartesian2D("test_PCG_P_small"))
1509 # suite.addTest(unittest.makeSuite(Test_StokesProblemCartesian3D))
1510 suite.addTest(unittest.makeSuite(Test_Rheologies))
1511 s=unittest.TextTestRunner(verbosity=2).run(suite)
1512 if not s.wasSuccessful(): sys.exit(1)
1513

  ViewVC Help
Powered by ViewVC 1.1.26