/[escript]/trunk/esys2/escript/py_src/linearPDEs.py
ViewVC logotype

Contents of /trunk/esys2/escript/py_src/linearPDEs.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 122 - (show annotations)
Thu Jun 9 05:38:05 2005 UTC (14 years, 4 months ago) by jgs
File MIME type: text/x-python
File size: 50734 byte(s)
Merge of development branch back to main trunk on 2005-06-09

1 # $Id$
2
3 ## @file linearPDEs.py
4
5 """
6 Functions and classes for linear PDEs
7 """
8
9 import escript
10 import util
11 import numarray
12
13
14 def _CompTuple2(t1,t2):
15 """
16 Compare two tuples
17
18 \param t1 The first tuple
19 \param t2 The second tuple
20 """
21
22 dif=t1[0]+t1[1]-(t2[0]+t2[1])
23 if dif<0: return 1
24 elif dif>0: return -1
25 else: return 0
26
27 def ELMAN_RAMAGE(P):
28 return (P-1.).wherePositive()*0.5*(1.-1./(P+1.e-15))
29
30 def SIMPLIFIED_BROOK_HUGHES(P):
31 c=(P-3.).whereNegative()
32 return P/6.*c+1./2.*(1.-c)
33
34 def HALF(P):
35 return escript.Scalar(0.5,P.getFunctionSpace())
36
37 class PDECoefficient:
38 """
39 A class for PDE coefficients
40 """
41 # identifier for location of Data objects defining COEFFICIENTS
42 INTERIOR=0
43 BOUNDARY=1
44 CONTACT=2
45 CONTINUOUS=3
46 # identifier in the pattern of COEFFICIENTS:
47 # the pattern is a tuple of EQUATION,SOLUTION,DIM where DIM represents the spatial dimension, EQUATION the number of equations and SOLUTION the
48 # number of unknowns.
49 EQUATION=3
50 SOLUTION=4
51 DIM=5
52 # indicator for what is altered if the coefficient is altered:
53 OPERATOR=5
54 RIGHTHANDSIDE=6
55 BOTH=7
56 def __init__(self,where,pattern,altering):
57 """
58 Initialise a PDE Coefficient type
59 """
60 self.what=where
61 self.pattern=pattern
62 self.altering=altering
63 self.resetValue()
64
65 def resetValue(self):
66 """
67 resets coefficient value to default
68 """
69 self.value=escript.Data()
70
71 def getFunctionSpace(self,domain):
72 """
73 defines the FunctionSpace of the coefficient on the domain
74
75 @param domain:
76 """
77 if self.what==self.INTERIOR: return escript.Function(domain)
78 elif self.what==self.BOUNDARY: return escript.FunctionOnBoundary(domain)
79 elif self.what==self.CONTACT: return escript.FunctionOnContactZero(domain)
80 elif self.what==self.CONTINUOUS: return escript.ContinuousFunction(domain)
81
82 def getValue(self):
83 """
84 returns the value of the coefficient:
85 """
86 return self.value
87
88 def setValue(self,newValue):
89 """
90 set the value of the coefficient to new value
91 """
92 self.value=newValue
93
94 def isAlteringOperator(self):
95 """
96 return true if the operator of the PDE is changed when the coefficient is changed
97 """
98 if self.altering==self.OPERATOR or self.altering==self.BOTH:
99 return not None
100 else:
101 return None
102
103 def isAlteringRightHandSide(self):
104 """
105 return true if the right hand side of the PDE is changed when the coefficient is changed
106 """
107 if self.altering==self.RIGHTHANDSIDE or self.altering==self.BOTH:
108 return not None
109 else:
110 return None
111
112 def estimateNumEquationsAndNumSolutions(self,shape=(),dim=3):
113 """
114 tries to estimate the number of equations in a given tensor shape for a given spatial dimension dim
115
116 @param shape:
117 @param dim:
118 """
119 if len(shape)>0:
120 num=max(shape)+1
121 else:
122 num=1
123 search=[]
124 for u in range(num):
125 for e in range(num):
126 search.append((e,u))
127 search.sort(_CompTuple2)
128 for item in search:
129 s=self.buildShape(item[0],item[1],dim)
130 if len(s)==0 and len(shape)==0:
131 return (1,1)
132 else:
133 if s==shape: return item
134 return None
135
136 def buildShape(self,e=1,u=1,dim=3):
137 """
138 builds the required shape for a given number of equations e, number of unknowns u and spatial dimension dim
139
140 @param e:
141 @param u:
142 @param dim:
143 """
144 s=()
145 for i in self.pattern:
146 if i==self.EQUATION:
147 if e>1: s=s+(e,)
148 elif i==self.SOLUTION:
149 if u>1: s=s+(u,)
150 else:
151 s=s+(dim,)
152 return s
153
154 class LinearPDE:
155 """
156 Class to handle a linear PDE
157
158 class to define a linear PDE of the form
159
160 \f[
161 -(A_{ijkl}u_{k,l})_{,j} -(B_{ijk}u_k)_{,j} + C_{ikl}u_{k,l} +D_{ik}u_k = - (X_{ij})_{,j} + Y_i
162 \f]
163
164 with boundary conditons:
165
166 \f[
167 n_j*(A_{ijkl}u_{k,l}+B_{ijk}u_k)_{,j} + d_{ik}u_k = - n_j*X_{ij} + y_i
168 \f]
169
170 and contact conditions
171
172 \f[
173 n_j*(A_{ijkl}u_{k,l}+B_{ijk}u_k)_{,j} + d_contact_{ik}[u_k] = - n_j*X_{ij} + y_contact_i
174 \f]
175
176 and constraints:
177
178 \f[
179 u_i=r_i \quad \mathrm{where} \quad q_i>0
180 \f]
181
182 """
183 TOL=1.e-13
184 DEFAULT_METHOD=util.DEFAULT_METHOD
185 DIRECT=util.DIRECT
186 CHOLEVSKY=util.CHOLEVSKY
187 PCG=util.PCG
188 CR=util.CR
189 CGS=util.CGS
190 BICGSTAB=util.BICGSTAB
191 SSOR=util.SSOR
192 GMRES=util.GMRES
193 PRES20=util.PRES20
194 ILU0=util.ILU0
195 JACOBI=util.JACOBI
196
197 def __init__(self,domain,numEquations=0,numSolutions=0):
198 """
199 initializes a new linear PDE.
200
201 @param args:
202 """
203 # COEFFICIENTS can be overwritten by subclasses:
204 self.COEFFICIENTS={
205 "A" : PDECoefficient(PDECoefficient.INTERIOR,(PDECoefficient.EQUATION,PDECoefficient.DIM,PDECoefficient.SOLUTION,PDECoefficient.DIM),PDECoefficient.OPERATOR),
206 "B" : PDECoefficient(PDECoefficient.INTERIOR,(PDECoefficient.EQUATION,PDECoefficient.DIM,PDECoefficient.SOLUTION),PDECoefficient.OPERATOR),
207 "C" : PDECoefficient(PDECoefficient.INTERIOR,(PDECoefficient.EQUATION,PDECoefficient.SOLUTION,PDECoefficient.DIM),PDECoefficient.OPERATOR),
208 "D" : PDECoefficient(PDECoefficient.INTERIOR,(PDECoefficient.EQUATION,PDECoefficient.SOLUTION),PDECoefficient.OPERATOR),
209 "X" : PDECoefficient(PDECoefficient.INTERIOR,(PDECoefficient.EQUATION,PDECoefficient.DIM),PDECoefficient.RIGHTHANDSIDE),
210 "Y" : PDECoefficient(PDECoefficient.INTERIOR,(PDECoefficient.EQUATION,),PDECoefficient.RIGHTHANDSIDE),
211 "d" : PDECoefficient(PDECoefficient.BOUNDARY,(PDECoefficient.EQUATION,PDECoefficient.SOLUTION),PDECoefficient.OPERATOR),
212 "y" : PDECoefficient(PDECoefficient.BOUNDARY,(PDECoefficient.EQUATION,),PDECoefficient.RIGHTHANDSIDE),
213 "d_contact" : PDECoefficient(PDECoefficient.CONTACT,(PDECoefficient.EQUATION,PDECoefficient.SOLUTION),PDECoefficient.OPERATOR),
214 "y_contact" : PDECoefficient(PDECoefficient.CONTACT,(PDECoefficient.EQUATION,),PDECoefficient.RIGHTHANDSIDE),
215 "r" : PDECoefficient(PDECoefficient.CONTINUOUS,(PDECoefficient.EQUATION,),PDECoefficient.RIGHTHANDSIDE),
216 "q" : PDECoefficient(PDECoefficient.CONTINUOUS,(PDECoefficient.SOLUTION,),PDECoefficient.BOTH)}
217
218 # initialize attributes
219 self.__debug=None
220 self.__domain=domain
221 self.__numEquations=numEquations
222 self.__numSolutions=numSolutions
223 self.cleanCoefficients()
224
225 self.__operator=escript.Operator()
226 self.__operator_isValid=False
227 self.__righthandside=escript.Data()
228 self.__righthandside_isValid=False
229 self.__solution=escript.Data()
230 self.__solution_isValid=False
231
232 # set some default values:
233 self.__homogeneous_constraint=True
234 self.__row_function_space=escript.Solution(self.__domain)
235 self.__column_function_space=escript.Solution(self.__domain)
236 self.__tolerance=1.e-8
237 self.__solver_method=util.DEFAULT_METHOD
238 self.__matrix_type=self.__domain.getSystemMatrixTypeId(util.DEFAULT_METHOD,False)
239 self.__sym=False
240 self.__lumping=False
241
242 def createCoefficient(self, name):
243 """
244 create a data object corresponding to coefficient name
245 @param name:
246 """
247 return escript.Data(shape = getShapeOfCoefficient(name), \
248 what = getFunctionSpaceForCoefficient(name))
249
250 def __del__(self):
251 pass
252
253 def getCoefficient(self,name):
254 """
255 return the value of the parameter name
256
257 @param name:
258 """
259 return self.COEFFICIENTS[name].getValue()
260
261 def getCoefficientOfPDE(self,name):
262 """
263 return the value of the coefficient name of the general PDE.
264 This method is called by the assembling routine it can be
265 overwritten to map coefficients of a particualr PDE to the general PDE.
266
267 @param name:
268 """
269 return self.getCoefficient(name)
270
271 def hasCoefficient(self,name):
272 """
273 return true if name is the name of a coefficient
274
275 @param name:
276 """
277 return self.COEFFICIENTS.has_key(name)
278
279 def getFunctionSpaceForEquation(self):
280 """
281 return true if the test functions should use reduced order
282 """
283 return self.__row_function_space
284
285 def getFunctionSpaceForSolution(self):
286 """
287 return true if the interpolation of the solution should use reduced order
288 """
289 return self.__column_function_space
290
291 def setValue(self,**coefficients):
292 """
293 sets new values to coefficients
294
295 @param coefficients:
296 """
297 self.__setValue(**coefficients)
298
299
300 def cleanCoefficients(self):
301 """
302 resets all coefficients to default values.
303 """
304 for i in self.COEFFICIENTS.iterkeys():
305 self.COEFFICIENTS[i].resetValue()
306
307 def createNewCoefficient(self,name):
308 """
309 returns a new coefficient appropriate for coefficient name:
310 """
311 return escript.Data(0,self.getShapeOfCoefficient(name),self.getFunctionSpaceForCoefficient(name))
312
313
314 def getShapeOfCoefficient(self,name):
315 """
316 return the shape of the coefficient name
317
318 @param name:
319 """
320 if self.hasCoefficient(name):
321 return self.COEFFICIENTS[name].buildShape(self.getNumEquations(),self.getNumSolutions(),self.getDomain().getDim())
322 else:
323 raise ValueError,"Solution coefficient %s requested"%name
324
325 def getFunctionSpaceForCoefficient(self,name):
326 """
327 return the atoms of the coefficient name
328
329 @param name:
330 """
331 if self.hasCoefficient(name):
332 return self.COEFFICIENTS[name].getFunctionSpace(self.getDomain())
333 else:
334 raise ValueError,"Solution coefficient %s requested"%name
335
336 def alteredCoefficient(self,name):
337 """
338 announce that coefficient name has been changed
339
340 @param name:
341 """
342 if self.hasCoefficient(name):
343 if self.COEFFICIENTS[name].isAlteringOperator(): self.__rebuildOperator()
344 if self.COEFFICIENTS[name].isAlteringRightHandSide(): self.__rebuildRightHandSide()
345 else:
346 raise ValueError,"unknown coefficient %s requested"%name
347
348 # ===== debug ==============================================================
349 def setDebugOn(self):
350 """
351 """
352 self.__debug=not None
353
354 def setDebugOff(self):
355 """
356 """
357 self.__debug=None
358
359 def debug(self):
360 """
361 returns true if the PDE is in the debug mode
362 """
363 return self.__debug
364
365 #===== Lumping ===========================
366 def setLumpingOn(self):
367 """
368 indicates to use matrix lumping
369 """
370 if not self.isUsingLumping():
371 if self.debug() : print "PDE Debug: lumping is set on"
372 self.__rebuildOperator()
373 self.__lumping=True
374
375 def setLumpingOff(self):
376 """
377 switches off matrix lumping
378 """
379 if self.isUsingLumping():
380 if self.debug() : print "PDE Debug: lumping is set off"
381 self.__rebuildOperator()
382 self.__lumping=False
383
384 def setLumping(self,flag=False):
385 """
386 set the matrix lumping flag to flag
387 """
388 if flag:
389 self.setLumpingOn()
390 else:
391 self.setLumpingOff()
392
393 def isUsingLumping(self):
394 """
395
396 """
397 return self.__lumping
398
399 #============ method business =========================================================
400 def setSolverMethod(self,solver=util.DEFAULT_METHOD):
401 """
402 sets a new solver
403 """
404 if not solver==self.getSolverMethod():
405 self.__solver_method=solver
406 if self.debug() : print "PDE Debug: New solver is %s"%solver
407 self.__checkMatrixType()
408
409 def getSolverMethod(self):
410 """
411 returns the solver method
412 """
413 return self.__solver_method
414
415 #============ tolerance business =========================================================
416 def setTolerance(self,tol=1.e-8):
417 """
418 resets the tolerance to tol.
419 """
420 if not tol>0:
421 raise ValueException,"Tolerance as to be positive"
422 if tol<self.getTolerance(): self.__rebuildSolution()
423 if self.debug() : print "PDE Debug: New tolerance %e",tol
424 self.__tolerance=tol
425 return
426 def getTolerance(self):
427 """
428 returns the tolerance set for the solution
429 """
430 return self.__tolerance
431
432 #===== symmetry flag ==========================
433 def isSymmetric(self):
434 """
435 returns true is the operator is considered to be symmetric
436 """
437 return self.__sym
438
439 def setSymmetryOn(self):
440 """
441 sets the symmetry flag to true
442 """
443 if not self.isSymmetric():
444 if self.debug() : print "PDE Debug: Operator is set to be symmetric"
445 self.__sym=True
446 self.__checkMatrixType()
447
448 def setSymmetryOff(self):
449 """
450 sets the symmetry flag to false
451 """
452 if self.isSymmetric():
453 if self.debug() : print "PDE Debug: Operator is set to be unsymmetric"
454 self.__sym=False
455 self.__checkMatrixType()
456
457 def setSymmetryTo(self,flag=False):
458 """
459 sets the symmetry flag to flag
460
461 @param flag:
462 """
463 if flag:
464 self.setSymmetryOn()
465 else:
466 self.setSymmetryOff()
467
468 #===== order reduction ==========================
469 def setReducedOrderOn(self):
470 """
471 switches to on reduced order
472 """
473 self.setReducedOrderForSolutionOn()
474 self.setReducedOrderForEquationOn()
475
476 def setReducedOrderOff(self):
477 """
478 switches to full order
479 """
480 self.setReducedOrderForSolutionOff()
481 self.setReducedOrderForEquationOff()
482
483 def setReducedOrderTo(self,flag=False):
484 """
485 sets order according to flag
486
487 @param flag:
488 """
489 self.setReducedOrderForSolutionTo(flag)
490 self.setReducedOrderForEquationTo(flag)
491
492
493 #===== order reduction solution ==========================
494 def setReducedOrderForSolutionOn(self):
495 """
496 switches to reduced order to interpolate solution
497 """
498 new_fs=escript.ReducedSolution(self.getDomain())
499 if self.getFunctionSpaceForSolution()!=new_fs:
500 if self.debug() : print "PDE Debug: Reduced order is used to interpolate solution."
501 self.__column_function_space=new_fs
502 self.__rebuildSystem(deep=True)
503
504 def setReducedOrderForSolutionOff(self):
505 """
506 switches to full order to interpolate solution
507 """
508 new_fs=escript.Solution(self.getDomain())
509 if self.getFunctionSpaceForSolution()!=new_fs:
510 if self.debug() : print "PDE Debug: Full order is used to interpolate solution."
511 self.__column_function_space=new_fs
512 self.__rebuildSystem(deep=True)
513
514 def setReducedOrderForSolutionTo(self,flag=False):
515 """
516 sets order for test functions according to flag
517
518 @param flag:
519 """
520 if flag:
521 self.setReducedOrderForSolutionOn()
522 else:
523 self.setReducedOrderForSolutionOff()
524
525 #===== order reduction equation ==========================
526 def setReducedOrderForEquationOn(self):
527 """
528 switches to reduced order for test functions
529 """
530 new_fs=escript.ReducedSolution(self.getDomain())
531 if self.getFunctionSpaceForEquation()!=new_fs:
532 if self.debug() : print "PDE Debug: Reduced order is used for test functions."
533 self.__row_function_space=new_fs
534 self.__rebuildSystem(deep=True)
535
536 def setReducedOrderForEquationOff(self):
537 """
538 switches to full order for test functions
539 """
540 new_fs=escript.Solution(self.getDomain())
541 if self.getFunctionSpaceForEquation()!=new_fs:
542 if self.debug() : print "PDE Debug: Full order is used for test functions."
543 self.__row_function_space=new_fs
544 self.__rebuildSystem(deep=True)
545
546 def setReducedOrderForEquationTo(self,flag=False):
547 """
548 sets order for test functions according to flag
549
550 @param flag:
551 """
552 if flag:
553 self.setReducedOrderForEquationOn()
554 else:
555 self.setReducedOrderForEquationOff()
556
557 # ==== initialization =====================================================================
558 def __getNewOperator(self):
559 """
560 """
561 return self.getDomain().newOperator( \
562 self.getNumEquations(), \
563 self.getFunctionSpaceForEquation(), \
564 self.getNumSolutions(), \
565 self.getFunctionSpaceForSolution(), \
566 self.__matrix_type)
567
568 def __makeFreshRightHandSide(self):
569 """
570 """
571 if self.debug() : print "PDE Debug: New right hand side allocated"
572 if self.getNumEquations()>1:
573 self.__righthandside=escript.Data(0.,(self.getNumEquations(),),self.getFunctionSpaceForEquation(),True)
574 else:
575 self.__righthandside=escript.Data(0.,(),self.getFunctionSpaceForEquation(),True)
576 return self.__righthandside
577
578 def __getNewSolution(self):
579 """
580 """
581 if self.debug() : print "PDE Debug: New right hand side allocated"
582 if self.getNumSolutions()>1:
583 return escript.Data(0.,(self.getNumSolutions(),),self.getFunctionSpaceForSolution(),True)
584 else:
585 return escript.Data(0.,(),self.getFunctionSpaceForSolution(),True)
586
587 def __makeFreshOperator(self):
588 """
589 """
590 if self.__operator.isEmpty():
591 self.__operator=self.__getNewOperator()
592 if self.debug() : print "PDE Debug: New operator allocated"
593 else:
594 self.__operator.setValue(0.)
595 self.__operator.resetSolver()
596 if self.debug() : print "PDE Debug: Operator reset to zero"
597 return self.__operator
598
599 #============ some serivice functions =====================================================
600 def getDomain(self):
601 """
602 returns the domain of the PDE
603 """
604 return self.__domain
605
606 def getDim(self):
607 """
608 returns the spatial dimension of the PDE
609 """
610 return self.getDomain().getDim()
611
612 def getNumEquations(self):
613 """
614 returns the number of equations
615 """
616 if self.__numEquations>0:
617 return self.__numEquations
618 else:
619 raise ValueError,"Number of equations is undefined. Please specify argument numEquations."
620
621 def getNumSolutions(self):
622 """
623 returns the number of unknowns
624 """
625 if self.__numSolutions>0:
626 return self.__numSolutions
627 else:
628 raise ValueError,"Number of solution is undefined. Please specify argument numSolutions."
629
630
631 def checkSymmetry(self,verbose=True):
632 """
633 returns if the Operator is symmetric. This is a very expensive operation!!! The symmetry flag is not altered.
634 """
635 verbose=verbose or self.debug()
636 out=True
637 if self.getNumSolutions()!=self.getNumEquations():
638 if verbose: print "non-symmetric PDE because of different number of equations and solutions"
639 out=False
640 else:
641 A=self.getCoefficientOfPDE("A")
642 if not A.isEmpty():
643 tol=util.Lsup(A)*self.TOL
644 if self.getNumSolutions()>1:
645 for i in range(self.getNumEquations()):
646 for j in range(self.getDim()):
647 for k in range(self.getNumSolutions()):
648 for l in range(self.getDim()):
649 if util.Lsup(A[i,j,k,l]-A[k,l,i,j])>tol:
650 if verbose: print "non-symmetric PDE because A[%d,%d,%d,%d]!=A[%d,%d,%d,%d]"%(i,j,k,l,k,l,i,j)
651 out=False
652 else:
653 for j in range(self.getDim()):
654 for l in range(self.getDim()):
655 if util.Lsup(A[j,l]-A[l,j])>tol:
656 if verbose: print "non-symmetric PDE because A[%d,%d]!=A[%d,%d]"%(j,l,l,j)
657 out=False
658 B=self.getCoefficientOfPDE("B")
659 C=self.getCoefficientOfPDE("C")
660 if B.isEmpty() and not C.isEmpty():
661 if verbose: print "non-symmetric PDE because B is not present but C is"
662 out=False
663 elif not B.isEmpty() and C.isEmpty():
664 if verbose: print "non-symmetric PDE because C is not present but B is"
665 out=False
666 elif not B.isEmpty() and not C.isEmpty():
667 tol=(util.Lsup(B)+util.Lsup(C))*self.TOL/2.
668 if self.getNumSolutions()>1:
669 for i in range(self.getNumEquations()):
670 for j in range(self.getDim()):
671 for k in range(self.getNumSolutions()):
672 if util.Lsup(B[i,j,k]-C[k,i,j])>tol:
673 if verbose: print "non-symmetric PDE because B[%d,%d,%d]!=C[%d,%d,%d]"%(i,j,k,k,i,j)
674 out=False
675 else:
676 for j in range(self.getDim()):
677 if util.Lsup(B[j]-C[j])>tol:
678 if verbose: print "non-symmetric PDE because B[%d]!=C[%d]"%(j,j)
679 out=False
680 if self.getNumSolutions()>1:
681 D=self.getCoefficientOfPDE("D")
682 if not D.isEmpty():
683 tol=util.Lsup(D)*self.TOL
684 for i in range(self.getNumEquations()):
685 for k in range(self.getNumSolutions()):
686 if util.Lsup(D[i,k]-D[k,i])>tol:
687 if verbose: print "non-symmetric PDE because D[%d,%d]!=D[%d,%d]"%(i,k,k,i)
688 out=False
689
690 return out
691
692 def getFlux(self,u):
693 """
694 returns the flux J_ij for a given u
695
696 \f[
697 J_ij=A_{ijkl}u_{k,l}+B_{ijk}u_k-X_{ij}
698 \f]
699
700 @param u: argument of the operator
701 """
702 raise SystemError,"getFlux is not implemented yet"
703 return None
704
705 def applyOperator(self,u):
706 """
707 applies the operator of the PDE to a given solution u in weak from
708
709 @param u: argument of the operator
710 """
711 return self.getOperator()*escript.Data(u,self.getFunctionSpaceForSolution())
712
713 def getResidual(self,u):
714 """
715 return the residual of u in the weak from
716
717 @param u:
718 """
719 return self.applyOperator(u)-self.getRightHandSide()
720
721 def __setValue(self,**coefficients):
722 """
723 sets new values to coefficient
724
725 @param coefficients:
726 """
727 # check if the coefficients are legal:
728 for i in coefficients.iterkeys():
729 if not self.hasCoefficient(i):
730 raise ValueError,"Attempt to set unknown coefficient %s"%i
731 # if the number of unknowns or equations is still unknown we try to estimate them:
732 if self.__numEquations<1 or self.__numSolutions<1:
733 for i,d in coefficients.iteritems():
734 if hasattr(d,"shape"):
735 s=d.shape
736 elif hasattr(d,"getShape"):
737 s=d.getShape()
738 else:
739 s=numarray.array(d).shape
740 if s!=None:
741 # get number of equations and number of unknowns:
742 res=self.COEFFICIENTS[i].estimateNumEquationsAndNumSolutions(s,self.getDim())
743 if res==None:
744 raise ValueError,"Illegal shape %s of coefficient %s"%(s,i)
745 else:
746 if self.__numEquations<1: self.__numEquations=res[0]
747 if self.__numSolutions<1: self.__numSolutions=res[1]
748 if self.__numEquations<1: raise ValueError,"unidententified number of equations"
749 if self.__numSolutions<1: raise ValueError,"unidententified number of solutions"
750 # now we check the shape of the coefficient if numEquations and numSolutions are set:
751 for i,d in coefficients.iteritems():
752 if d==None:
753 d2=escript.Data()
754 elif isinstance(d,escript.Data):
755 if d.isEmpty():
756 d2=d
757 else:
758 d2=escript.Data(d,self.getFunctionSpaceForCoefficient(i))
759 else:
760 d2=escript.Data(d,self.getFunctionSpaceForCoefficient(i))
761 if not d2.isEmpty():
762 if not self.getShapeOfCoefficient(i)==d2.getShape():
763 raise ValueError,"Expected shape for coefficient %s is %s but actual shape is %s."%(i,self.getShapeOfCoefficient(i),d2.getShape())
764 # overwrite new values:
765 if self.debug(): print "PDE Debug: Coefficient %s has been altered."%i
766 self.COEFFICIENTS[i].setValue(d2)
767 self.alteredCoefficient(i)
768
769 # reset the HomogeneousConstraintFlag:
770 self.__setHomogeneousConstraintFlag()
771 if len(coefficients)>0 and not self.isUsingLumping() and not self.__homogeneous_constraint: self.__rebuildSystem()
772
773 def __setHomogeneousConstraintFlag(self):
774 """
775 checks if the constraints are homogeneous and sets self.__homogeneous_constraint accordingly.
776 """
777 self.__homogeneous_constraint=True
778 q=self.getCoefficientOfPDE("q")
779 r=self.getCoefficientOfPDE("r")
780 if not q.isEmpty() and not r.isEmpty():
781 if (q*r).Lsup()>=1.e-13*r.Lsup(): self.__homogeneous_constraint=False
782 if self.debug():
783 if self.__homogeneous_constraint:
784 print "PDE Debug: Constraints are homogeneous."
785 else:
786 print "PDE Debug: Constraints are inhomogeneous."
787
788
789 # ==== rebuild switches =====================================================================
790 def __rebuildSolution(self,deep=False):
791 """
792 indicates the PDE has to be reolved if the solution is requested
793 """
794 if self.__solution_isValid and self.debug() : print "PDE Debug: PDE has to be resolved."
795 self.__solution_isValid=False
796 if deep: self.__solution=escript.Data()
797
798
799 def __rebuildOperator(self,deep=False):
800 """
801 indicates the operator has to be rebuilt next time it is used
802 """
803 if self.__operator_isValid and self.debug() : print "PDE Debug: Operator has to be rebuilt."
804 self.__rebuildSolution(deep)
805 self.__operator_isValid=False
806 if deep: self.__operator=escript.Operator()
807
808 def __rebuildRightHandSide(self,deep=False):
809 """
810 indicates the right hand side has to be rebuild next time it is used
811 """
812 if self.__righthandside_isValid and self.debug() : print "PDE Debug: Right hand side has to be rebuilt."
813 self.__rebuildSolution(deep)
814 self.__righthandside_isValid=False
815 if deep: self.__righthandside=escript.Data()
816
817 def __rebuildSystem(self,deep=False):
818 """
819 annonced that all coefficient name has been changed
820 """
821 self.__rebuildSolution(deep)
822 self.__rebuildOperator(deep)
823 self.__rebuildRightHandSide(deep)
824
825 def __checkMatrixType(self):
826 """
827 reassess the matrix type and, if needed, initiates an operator rebuild
828 """
829 new_matrix_type=self.getDomain().getSystemMatrixTypeId(self.getSolverMethod(),self.isSymmetric())
830 if not new_matrix_type==self.__matrix_type:
831 if self.debug() : print "PDE Debug: Matrix type is now %d."%new_matrix_type
832 self.__matrix_type=new_matrix_type
833 self.__rebuildOperator(deep=True)
834
835 #============ assembling =======================================================
836 def __copyConstraint(self):
837 """
838 copies the constrint condition into u
839 """
840 if not self.__righthandside.isEmpty():
841 q=self.getCoefficientOfPDE("q")
842 r=self.getCoefficientOfPDE("r")
843 if not q.isEmpty():
844 if r.isEmpty():
845 r2=escript.Data(0,self.__righthandside.getShape(),self.__righthandside.getFunctionSpace())
846 else:
847 r2=escript.Data(r,self.__righthandside.getFunctionSpace())
848 self.__righthandside.copyWithMask(r2,escript.Data(q,self.__righthandside.getFunctionSpace()))
849
850 def __applyConstraint(self):
851 """
852 applies the constraints defined by q and r to the system
853 """
854 q=self.getCoefficientOfPDE("q")
855 r=self.getCoefficientOfPDE("r")
856 if not q.isEmpty() and not self.__operator.isEmpty():
857 # q is the row and column mask to indicate where constraints are set:
858 row_q=escript.Data(q,self.getFunctionSpaceForEquation())
859 col_q=escript.Data(q,self.getFunctionSpaceForSolution())
860 u=self.__getNewSolution()
861 if r.isEmpty():
862 r_s=self.__getNewSolution()
863 else:
864 r_s=escript.Data(r,self.getFunctionSpaceForSolution())
865 u.copyWithMask(r_s,col_q)
866 if self.isUsingLumping():
867 self.__operator.copyWithMask(escript.Data(1,q.getShape(),self.getFunctionSpaceForEquation()),row_q)
868 else:
869 if not self.__righthandside.isEmpty(): self.__righthandside-=self.__operator*u
870 self.__operator.nullifyRowsAndCols(row_q,col_q,1.)
871
872 def getSystem(self):
873 """
874 return the operator and right hand side of the PDE
875 """
876 if not self.__operator_isValid or not self.__righthandside_isValid:
877 if self.isUsingLumping():
878 if not self.__operator_isValid:
879 if not self.getFunctionSpaceForEquation()==self.getFunctionSpaceForSolution():
880 raise TypeError,"Lumped matrix requires same order for equations and unknowns"
881 if not self.getCoefficientOfPDE("A").isEmpty():
882 raise Warning,"Lumped matrix does not allow coefficient A"
883 if not self.getCoefficientOfPDE("B").isEmpty():
884 raise Warning,"Lumped matrix does not allow coefficient B"
885 if not self.getCoefficientOfPDE("C").isEmpty():
886 raise Warning,"Lumped matrix does not allow coefficient C"
887 if self.debug() : print "PDE Debug: New lumped operator is built."
888 mat=self.__getNewOperator()
889 self.getDomain().addPDEToSystem(mat,escript.Data(), \
890 self.getCoefficientOfPDE("A"), \
891 self.getCoefficientOfPDE("B"), \
892 self.getCoefficientOfPDE("C"), \
893 self.getCoefficientOfPDE("D"), \
894 escript.Data(), \
895 escript.Data(), \
896 self.getCoefficientOfPDE("d"), \
897 escript.Data(),\
898 self.getCoefficientOfPDE("d_contact"), \
899 escript.Data())
900 self.__operator=mat*escript.Data(1,(self.getNumSolutions(),),self.getFunctionSpaceForSolution(),True)
901 self.__applyConstraint()
902 self.__operator_isValid=True
903 if not self.__righthandside_isValid:
904 if self.debug() : print "PDE Debug: New right hand side is built."
905 self.getDomain().addPDEToRHS(self.__makeFreshRightHandSide(), \
906 self.getCoefficientOfPDE("X"), \
907 self.getCoefficientOfPDE("Y"),\
908 self.getCoefficientOfPDE("y"),\
909 self.getCoefficientOfPDE("y_contact"))
910 self.__copyConstraint()
911 self.__righthandside_isValid=True
912 else:
913 if not self.__operator_isValid and not self.__righthandside_isValid:
914 if self.debug() : print "PDE Debug: New system is built."
915 self.getDomain().addPDEToSystem(self.__makeFreshOperator(),self.__makeFreshRightHandSide(), \
916 self.getCoefficientOfPDE("A"), \
917 self.getCoefficientOfPDE("B"), \
918 self.getCoefficientOfPDE("C"), \
919 self.getCoefficientOfPDE("D"), \
920 self.getCoefficientOfPDE("X"), \
921 self.getCoefficientOfPDE("Y"), \
922 self.getCoefficientOfPDE("d"), \
923 self.getCoefficientOfPDE("y"), \
924 self.getCoefficientOfPDE("d_contact"), \
925 self.getCoefficientOfPDE("y_contact"))
926 self.__applyConstraint()
927 self.__copyConstraint()
928 self.__operator_isValid=True
929 self.__righthandside_isValid=True
930 elif not self.__righthandside_isValid:
931 if self.debug() : print "PDE Debug: New right hand side is built."
932 self.getDomain().addPDEToRHS(self.__makeFreshRightHandSide(), \
933 self.getCoefficientOfPDE("X"), \
934 self.getCoefficientOfPDE("Y"),\
935 self.getCoefficientOfPDE("y"),\
936 self.getCoefficientOfPDE("y_contact"))
937 self.__copyConstraint()
938 self.__righthandside_isValid=True
939 elif not self.__operator_isValid:
940 if self.debug() : print "PDE Debug: New operator is built."
941 self.getDomain().addPDEToSystem(self.__makeFreshOperator(),escript.Data(), \
942 self.getCoefficientOfPDE("A"), \
943 self.getCoefficientOfPDE("B"), \
944 self.getCoefficientOfPDE("C"), \
945 self.getCoefficientOfPDE("D"), \
946 escript.Data(), \
947 escript.Data(), \
948 self.getCoefficientOfPDE("d"), \
949 escript.Data(),\
950 self.getCoefficientOfPDE("d_contact"), \
951 escript.Data())
952 self.__applyConstraint()
953 self.__operator_isValid=True
954 return (self.__operator,self.__righthandside)
955 def getOperator(self):
956 """
957 returns the operator of the PDE
958 """
959 return self.getSystem()[0]
960
961 def getRightHandSide(self):
962 """
963 returns the right hand side of the PDE
964 """
965 return self.getSystem()[1]
966
967 def solve(self,**options):
968 """
969 solve the PDE
970
971 @param options:
972 """
973 mat,f=self.getSystem()
974 if self.isUsingLumping():
975 out=f/mat
976 else:
977 options[util.TOLERANCE_KEY]=self.getTolerance()
978 options[util.METHOD_KEY]=self.getSolverMethod()
979 options[util.SYMMETRY_KEY]=self.isSymmetric()
980 if self.debug() : print "PDE Debug: solver options: ",options
981 out=mat.solve(f,options)
982 return out
983
984 def getSolution(self,**options):
985 """
986 returns the solution of the PDE
987
988 @param options:
989 """
990 if not self.__solution_isValid:
991 if self.debug() : print "PDE Debug: PDE is resolved."
992 self.__solution=self.solve(**options)
993 self.__solution_isValid=True
994 return self.__solution
995
996
997
998 def ELMAN_RAMAGE(P):
999 """ """
1000 return (P-1.).wherePositive()*0.5*(1.-1./(P+1.e-15))
1001 def SIMPLIFIED_BROOK_HUGHES(P):
1002 """ """
1003 c=(P-3.).whereNegative()
1004 return P/6.*c+1./2.*(1.-c)
1005
1006 def HALF(P):
1007 """ """
1008 return escript.Scalar(0.5,P.getFunctionSpace())
1009
1010 class AdvectivePDE(LinearPDE):
1011 """
1012 Class to handle a linear PDE dominated by advective terms:
1013
1014 class to define a linear PDE of the form
1015
1016 \f[
1017 -(A_{ijkl}u_{k,l})_{,j} -(B_{ijk}u_k)_{,j} + C_{ikl}u_{k,l} +D_{ik}u_k = - (X_{ij})_{,j} + Y_i
1018 \f]
1019
1020 with boundary conditons:
1021
1022 \f[
1023 n_j*(A_{ijkl}u_{k,l}+B_{ijk}u_k)_{,j} + d_{ik}u_k = - n_j*X_{ij} + y_i
1024 \f]
1025
1026 and contact conditions
1027
1028 \f[
1029 n_j*(A_{ijkl}u_{k,l}+B_{ijk}u_k)_{,j} + d^{contact}_{ik}[u_k] = - n_j*X_{ij} + y^{contact}_{i}
1030 \f]
1031
1032 and constraints:
1033
1034 \f[
1035 u_i=r_i \quad \mathrm{where} \quad q_i>0
1036 \f]
1037 """
1038 def __init__(self,domain,numEquations=0,numSolutions=0,xi=ELMAN_RAMAGE):
1039 LinearPDE.__init__(self,domain,numEquations,numSolutions)
1040 self.__xi=xi
1041 self.__Xi=escript.Data()
1042
1043 def __calculateXi(self,peclet_factor,Z,h):
1044 Z_max=util.Lsup(Z)
1045 if Z_max>0.:
1046 return h*self.__xi(Z*peclet_factor)/(Z+Z_max*self.TOL)
1047 else:
1048 return 0.
1049
1050 def setValue(self,**args):
1051 if "A" in args.keys() or "B" in args.keys() or "C" in args.keys(): self.__Xi=escript.Data()
1052 self._LinearPDE__setValue(**args)
1053
1054 def getXi(self):
1055 if self.__Xi.isEmpty():
1056 B=self.getCoefficient("B")
1057 C=self.getCoefficient("C")
1058 A=self.getCoefficient("A")
1059 h=self.getDomain().getSize()
1060 self.__Xi=escript.Scalar(0.,self.getFunctionSpaceForCoefficient("A"))
1061 if not C.isEmpty() or not B.isEmpty():
1062 if not C.isEmpty() and not B.isEmpty():
1063 Z2=escript.Scalar(0,self.getFunctionSpaceForCoefficient("A"))
1064 if self.getNumEquations()>1:
1065 if self.getNumSolutions()>1:
1066 for i in range(self.getNumEquations()):
1067 for k in range(self.getNumSolutions()):
1068 for l in range(self.getDim()): Z2+=(C[i,k,l]-B[i,l,k])**2
1069 else:
1070 for i in range(self.getNumEquations()):
1071 for l in range(self.getDim()): Z2+=(C[i,l]-B[i,l])**2
1072 else:
1073 if self.getNumSolutions()>1:
1074 for k in range(self.getNumSolutions()):
1075 for l in range(self.getDim()): Z2+=(C[k,l]-B[l,k])**2
1076 else:
1077 for l in range(self.getDim()): Z2+=(C[l]-B[l])**2
1078 length_of_Z=util.sqrt(Z2)
1079 elif C.isEmpty():
1080 length_of_Z=util.length(B)
1081 else:
1082 length_of_Z=util.length(C)
1083
1084 Z_max=util.Lsup(length_of_Z)
1085 if Z_max>0.:
1086 length_of_A=util.length(A)
1087 A_max=util.Lsup(length_of_A)
1088 if A_max>0:
1089 inv_A=1./(length_of_A+A_max*self.TOL)
1090 else:
1091 inv_A=1./self.TOL
1092 peclet_number=length_of_Z*h/2*inv_A
1093 xi=self.__xi(peclet_number)
1094 self.__Xi=h*xi/(length_of_Z+Z_max*self.TOL)
1095 print "@ preclet number = %e"%util.Lsup(peclet_number),util.Lsup(xi),util.Lsup(length_of_Z)
1096 return self.__Xi
1097
1098
1099 def getCoefficientOfPDE(self,name):
1100 """
1101 return the value of the coefficient name of the general PDE
1102
1103 @param name:
1104 """
1105 if not self.getNumEquations() == self.getNumSolutions():
1106 raise ValueError,"AdvectivePDE expects the number of solution componets and the number of equations to be equal."
1107
1108 if name == "A" :
1109 A=self.getCoefficient("A")
1110 B=self.getCoefficient("B")
1111 C=self.getCoefficient("C")
1112 if B.isEmpty() and C.isEmpty():
1113 Aout=A
1114 else:
1115 if A.isEmpty():
1116 Aout=self.createNewCoefficient("A")
1117 else:
1118 Aout=A[:]
1119 Xi=self.getXi()
1120 if self.getNumEquations()>1:
1121 for i in range(self.getNumEquations()):
1122 for j in range(self.getDim()):
1123 for k in range(self.getNumSolutions()):
1124 for l in range(self.getDim()):
1125 if not C.isEmpty() and not B.isEmpty():
1126 for p in range(self.getNumEquations()): Aout[i,j,k,l]+=Xi*(C[p,i,j]-B[p,j,i])*(C[p,k,l]-B[p,l,k])
1127 elif C.isEmpty():
1128 for p in range(self.getNumEquations()): Aout[i,j,k,l]+=Xi*B[p,j,i]*B[p,l,k]
1129 else:
1130 for p in range(self.getNumEquations()): Aout[i,j,k,l]+=Xi*C[p,i,j]*C[p,k,l]
1131 else:
1132 for j in range(self.getDim()):
1133 for l in range(self.getDim()):
1134 if not C.isEmpty() and not B.isEmpty():
1135 Aout[j,l]+=Xi*(C[j]-B[j])*(C[l]-B[l])
1136 elif C.isEmpty():
1137 Aout[j,l]+=Xi*B[j]*B[l]
1138 else:
1139 Aout[j,l]+=Xi*C[j]*C[l]
1140 return Aout
1141 elif name == "B" :
1142 B=self.getCoefficient("B")
1143 C=self.getCoefficient("C")
1144 D=self.getCoefficient("D")
1145 if C.isEmpty() or D.isEmpty():
1146 Bout=B
1147 else:
1148 Xi=self.getXi()
1149 if B.isEmpty():
1150 Bout=self.createNewCoefficient("B")
1151 else:
1152 Bout=B[:]
1153 if self.getNumEquations()>1:
1154 for k in range(self.getNumSolutions()):
1155 for p in range(self.getNumEquations()):
1156 tmp=Xi*D[p,k]
1157 for i in range(self.getNumEquations()):
1158 for j in range(self.getDim()):
1159 Bout[i,j,k]+=tmp*C[p,i,j]
1160 else:
1161 tmp=Xi*D
1162 for j in range(self.getDim()): Bout[j]+=tmp*C[j]
1163 return Bout
1164 elif name == "C" :
1165 B=self.getCoefficient("B")
1166 C=self.getCoefficient("C")
1167 D=self.getCoefficient("D")
1168 if B.isEmpty() or D.isEmpty():
1169 Cout=C
1170 else:
1171 Xi=self.getXi()
1172 if C.isEmpty():
1173 Cout=self.createNewCoefficient("C")
1174 else:
1175 Cout=C[:]
1176 if self.getNumEquations()>1:
1177 for k in range(self.getNumSolutions()):
1178 for p in range(self.getNumEquations()):
1179 tmp=Xi*D[p,k]
1180 for i in range(self.getNumEquations()):
1181 for l in range(self.getDim()):
1182 Cout[i,k,l]+=tmp*B[p,l,i]
1183 else:
1184 tmp=Xi*D
1185 for j in range(self.getDim()): Cout[j]+=tmp*B[j]
1186 return Cout
1187 elif name == "D" :
1188 return self.getCoefficient("D")
1189 elif name == "X" :
1190 X=self.getCoefficient("X")
1191 Y=self.getCoefficient("Y")
1192 B=self.getCoefficient("B")
1193 C=self.getCoefficient("C")
1194 if Y.isEmpty() or (B.isEmpty() and C.isEmpty()):
1195 Xout=X
1196 else:
1197 if X.isEmpty():
1198 Xout=self.createNewCoefficient("X")
1199 else:
1200 Xout=X[:]
1201 Xi=self.getXi()
1202 if self.getNumEquations()>1:
1203 for p in range(self.getNumEquations()):
1204 tmp=Xi*Y[p]
1205 for i in range(self.getNumEquations()):
1206 for j in range(self.getDim()):
1207 if not C.isEmpty() and not B.isEmpty():
1208 Xout[i,j]+=tmp*(C[p,i,j]-B[p,j,i])
1209 elif C.isEmpty():
1210 Xout[i,j]-=tmp*B[p,j,i]
1211 else:
1212 Xout[i,j]+=tmp*C[p,i,j]
1213 else:
1214 tmp=Xi*Y
1215 for j in range(self.getDim()):
1216 if not C.isEmpty() and not B.isEmpty():
1217 Xout[j]+=tmp*(C[j]-B[j])
1218 elif C.isEmpty():
1219 Xout[j]-=tmp*B[j]
1220 else:
1221 Xout[j]+=tmp*C[j]
1222 return Xout
1223 elif name == "Y" :
1224 return self.getCoefficient("Y")
1225 elif name == "d" :
1226 return self.getCoefficient("d")
1227 elif name == "y" :
1228 return self.getCoefficient("y")
1229 elif name == "d_contact" :
1230 return self.getCoefficient("d_contact")
1231 elif name == "y_contact" :
1232 return self.getCoefficient("y_contact")
1233 elif name == "r" :
1234 return self.getCoefficient("r")
1235 elif name == "q" :
1236 return self.getCoefficient("q")
1237 else:
1238 raise SystemError,"unknown PDE coefficient %s",name
1239
1240
1241 class Poisson(LinearPDE):
1242 """
1243 Class to define a Poisson equstion problem:
1244
1245 class to define a linear PDE of the form
1246 \f[
1247 -u_{,jj} = f
1248 \f]
1249
1250 with boundary conditons:
1251
1252 \f[
1253 n_j*u_{,j} = 0
1254 \f]
1255
1256 and constraints:
1257
1258 \f[
1259 u=0 \quad \mathrm{where} \quad q>0
1260 \f]
1261 """
1262
1263 def __init__(self,domain,f=escript.Data(),q=escript.Data()):
1264 LinearPDE.__init__(self,domain,1,1)
1265 self.COEFFICIENTS={
1266 "f" : PDECoefficient(PDECoefficient.INTERIOR,(PDECoefficient.EQUATION,),PDECoefficient.RIGHTHANDSIDE),
1267 "q" : PDECoefficient(PDECoefficient.CONTINUOUS,(PDECoefficient.EQUATION,),PDECoefficient.BOTH)}
1268 self.setSymmetryOn()
1269 self.setValue(f,q)
1270
1271 def setValue(self,f=escript.Data(),q=escript.Data()):
1272 self._LinearPDE__setValue(f=f,q=q)
1273
1274 def getCoefficientOfPDE(self,name):
1275 """
1276 return the value of the coefficient name of the general PDE
1277
1278 @param name:
1279 """
1280 if name == "A" :
1281 return escript.Data(numarray.identity(self.getDim()),escript.Function(self.getDomain()))
1282 elif name == "B" :
1283 return escript.Data()
1284 elif name == "C" :
1285 return escript.Data()
1286 elif name == "D" :
1287 return escript.Data()
1288 elif name == "X" :
1289 return escript.Data()
1290 elif name == "Y" :
1291 return self.getCoefficient("f")
1292 elif name == "d" :
1293 return escript.Data()
1294 elif name == "y" :
1295 return escript.Data()
1296 elif name == "d_contact" :
1297 return escript.Data()
1298 elif name == "y_contact" :
1299 return escript.Data()
1300 elif name == "r" :
1301 return escript.Data()
1302 elif name == "q" :
1303 return self.getCoefficient("q")
1304 else:
1305 raise SystemError,"unknown PDE coefficient %s",name
1306
1307 # $Log$
1308 # Revision 1.8 2005/06/09 05:37:59 jgs
1309 # Merge of development branch back to main trunk on 2005-06-09
1310 #
1311 # Revision 1.7 2005/05/06 04:26:10 jgs
1312 # Merge of development branch back to main trunk on 2005-05-06
1313 #
1314 # Revision 1.1.2.23 2005/05/13 00:55:20 cochrane
1315 # Fixed up some docstrings. Moved module-level functions to top of file so
1316 # that epydoc and doxygen can pick them up properly.
1317 #
1318 # Revision 1.1.2.22 2005/05/12 11:41:30 gross
1319 # some basic Models have been added
1320 #
1321 # Revision 1.1.2.21 2005/05/12 07:16:12 cochrane
1322 # Moved ELMAN_RAMAGE, SIMPLIFIED_BROOK_HUGHES, and HALF functions to bottom of
1323 # file so that the AdvectivePDE class is picked up by doxygen. Some
1324 # reformatting of docstrings. Addition of code to make equations come out
1325 # as proper LaTeX.
1326 #
1327 # Revision 1.1.2.20 2005/04/15 07:09:08 gross
1328 # some problems with functionspace and linearPDEs fixed.
1329 #
1330 # Revision 1.1.2.19 2005/03/04 05:27:07 gross
1331 # bug in SystemPattern fixed.
1332 #
1333 # Revision 1.1.2.18 2005/02/08 06:16:45 gross
1334 # Bugs in AdvectivePDE fixed, AdvectiveTest is stable but more testing is needed
1335 #
1336 # Revision 1.1.2.17 2005/02/08 05:56:19 gross
1337 # Reference Number handling added
1338 #
1339 # Revision 1.1.2.16 2005/02/07 04:41:28 gross
1340 # some function exposed to python to make mesh merging running
1341 #
1342 # Revision 1.1.2.15 2005/02/03 00:14:44 gross
1343 # timeseries add and ESySParameter.py renames esysXML.py for consistence
1344 #
1345 # Revision 1.1.2.14 2005/02/01 06:44:10 gross
1346 # new implementation of AdvectivePDE which now also updates right hand side. systems of PDEs are still not working
1347 #
1348 # Revision 1.1.2.13 2005/01/25 00:47:07 gross
1349 # updates in the documentation
1350 #
1351 # Revision 1.1.2.12 2005/01/12 01:28:04 matt
1352 # Added createCoefficient method for linearPDEs.
1353 #
1354 # Revision 1.1.2.11 2005/01/11 01:55:34 gross
1355 # a problem in linearPDE class fixed
1356 #
1357 # Revision 1.1.2.10 2005/01/07 01:13:29 gross
1358 # some bugs in linearPDE fixed
1359 #
1360 # Revision 1.1.2.9 2005/01/06 06:24:58 gross
1361 # some bugs in slicing fixed
1362 #
1363 # Revision 1.1.2.8 2005/01/05 04:21:40 gross
1364 # FunctionSpace checking/matchig in slicing added
1365 #
1366 # Revision 1.1.2.7 2004/12/29 10:03:41 gross
1367 # bug in setValue fixed
1368 #
1369 # Revision 1.1.2.6 2004/12/29 05:29:59 gross
1370 # AdvectivePDE successfully tested for Peclet number 1000000. there is still a problem with setValue and Data()
1371 #
1372 # Revision 1.1.2.5 2004/12/29 00:18:41 gross
1373 # AdvectivePDE added
1374 #
1375 # Revision 1.1.2.4 2004/12/24 06:05:41 gross
1376 # some changes in linearPDEs to add AdevectivePDE
1377 #
1378 # Revision 1.1.2.3 2004/12/16 00:12:34 gross
1379 # __init__ of LinearPDE does not accept any coefficient anymore
1380 #
1381 # Revision 1.1.2.2 2004/12/14 03:55:01 jgs
1382 # *** empty log message ***
1383 #
1384 # Revision 1.1.2.1 2004/12/12 22:53:47 gross
1385 # linearPDE has been renamed LinearPDE
1386 #
1387 # Revision 1.1.1.1.2.7 2004/12/07 10:13:08 gross
1388 # GMRES added
1389 #
1390 # Revision 1.1.1.1.2.6 2004/12/07 03:19:50 gross
1391 # options for GMRES and PRES20 added
1392 #
1393 # Revision 1.1.1.1.2.5 2004/12/01 06:25:15 gross
1394 # some small changes
1395 #
1396 # Revision 1.1.1.1.2.4 2004/11/24 01:50:21 gross
1397 # Finley solves 4M unknowns now
1398 #
1399 # Revision 1.1.1.1.2.3 2004/11/15 06:05:26 gross
1400 # poisson solver added
1401 #
1402 # Revision 1.1.1.1.2.2 2004/11/12 06:58:15 gross
1403 # a lot of changes to get the linearPDE class running: most important change is that there is no matrix format exposed to the user anymore. the format is chosen by the Domain according to the solver and symmetry
1404 #
1405 # Revision 1.1.1.1.2.1 2004/10/28 22:59:22 gross
1406 # finley's RecTest.py is running now: problem in SystemMatrixAdapater fixed
1407 #
1408 # Revision 1.1.1.1 2004/10/26 06:53:56 jgs
1409 # initial import of project esys2
1410 #
1411 # Revision 1.3.2.3 2004/10/26 06:43:48 jgs
1412 # committing Lutz's and Paul's changes to brach jgs
1413 #
1414 # Revision 1.3.4.1 2004/10/20 05:32:51 cochrane
1415 # Added incomplete Doxygen comments to files, or merely put the docstrings that already exist into Doxygen form.
1416 #
1417 # Revision 1.3 2004/09/23 00:53:23 jgs
1418 # minor fixes
1419 #
1420 # Revision 1.1 2004/08/28 12:58:06 gross
1421 # SimpleSolve is not running yet: problem with == of functionsspace
1422 #
1423 #

Properties

Name Value
svn:eol-style native
svn:keywords Author Date Id Revision

  ViewVC Help
Powered by ViewVC 1.1.26