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

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

Parent Directory Parent Directory | Revision Log Revision Log


Revision 121 - (show annotations)
Fri May 6 04:26:16 2005 UTC (14 years, 7 months ago) by jgs
Original Path: trunk/esys2/escript/py_src/linearPDEs.py
File MIME type: text/x-python
File size: 51196 byte(s)
Merge of development branch back to main trunk on 2005-05-06

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

Properties

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

  ViewVC Help
Powered by ViewVC 1.1.26