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

Annotation of /trunk/finley/test/python/run_models.py

Parent Directory Parent Directory | Revision Log Revision Log


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