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

  ViewVC Help
Powered by ViewVC 1.1.26