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

Annotation of /trunk/escript/py_src/rheologies.py

Parent Directory Parent Directory | Revision Log Revision Log


Revision 2300 - (hide annotations)
Wed Mar 11 08:17:57 2009 UTC (11 years, 5 months ago) by gross
File MIME type: text/x-python
File size: 22314 byte(s)
a power law model. tests still need to be added and needs to be linked with an incompressible solver.
1 ksteube 1809
2     ########################################################
3 gross 1639 #
4 ksteube 1809 # Copyright (c) 2003-2008 by University of Queensland
5     # Earth Systems Science Computational Center (ESSCC)
6     # http://www.uq.edu.au/esscc
7 gross 1639 #
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 1639 #
12 ksteube 1809 ########################################################
13 gross 1639
14 ksteube 1809 __copyright__="""Copyright (c) 2003-2008 by University of Queensland
15     Earth Systems Science Computational Center (ESSCC)
16     http://www.uq.edu.au/esscc
17     Primary Business: Queensland, Australia"""
18     __license__="""Licensed under the Open Software License version 3.0
19     http://www.opensource.org/licenses/osl-3.0.php"""
20     __url__="http://www.uq.edu.au/esscc/escript-finley"
21    
22 gross 1639 """
23     Some models for flow
24    
25     @var __author__: name of author
26     @var __copyright__: copyrights
27     @var __license__: licence agreement
28     @var __url__: url entry point on documentation
29     @var __version__: version
30     @var __date__: date of the version
31     """
32    
33     __author__="Lutz Gross, l.gross@uq.edu.au"
34    
35     from escript import *
36     import util
37     from linearPDEs import LinearPDE
38 gross 2100 from pdetools import Defect, NewtonGMRES, ArithmeticTuple
39 gross 1639
40 gross 2300 class PowerLaw(object):
41     """
42     this implements the power law for a composition of a set of materials where the viscosity eta of each material is given by a
43     power law relationship of the form
44    
45     M{eta=eta_N*(tau/tau_t)**(1.-1./power)}
46    
47     where tau is equivalent stress and eta_N, tau_t and power are given constant. Moreover an elastic component can be considered.
48     Moreover tau meets the Drucker-Prager type yield condition
49    
50     M{tau <= tau_Y + friction * pressure}
51    
52     where gamma_dot is the equivalent.
53     """
54     def __init__(self, numMaterials=1,verbose=False):
55     """
56     initializes a power law
57    
58     @param numMaterials: number of materials
59     @type numMaterials: C{int}
60     @param verbose: if C{True} some informations are printed.
61     @type verbose: C{bool}
62     """
63     if numMaterials<1:
64     raise ValueError,"at least one material must be defined."
65     self.__numMaterials=numMaterials
66     self.__eta_N=[None for i in xrange(self.__numMaterials)]
67     self.__tau_t=[1. for i in xrange(self.__numMaterials)]
68     self.__power=[1. for i in xrange(self.__numMaterials)]
69     self.__tau_Y=None
70     self.__friction=0
71     self.__mu=None
72     self.__verbose=verbose
73     self.setEtaTolerance()
74     #===========================================================================
75     def getNumMaterials(self):
76     """
77     returns the numebr of materials
78    
79     @return: number of materials
80     @rtype: C{int}
81     """
82     return self.__numMaterials
83     def validMaterialId(self,id=0):
84     """
85     checks if a given material id is valid
86    
87     @param id: a material id
88     @type id: C{int}
89     @return: C{True} is the id is valid
90     @rtype: C{bool}
91     """
92     return 0<=id and id<self.getNumMaterials()
93     def setEtaTolerance(self,rtol=util.sqrt(util.EPSILON)):
94     """
95     sets the relative tolerance for the effectice viscosity.
96    
97     @param rtol: relative tolerance
98     @type rtol: positive C{float}
99     """
100     if rtol<=0:
101     raise ValueError,"rtol needs to positive."
102     self.__rtol=rtol
103     def getEtaTolerance(self):
104     """
105     returns the relative tolerance for the effectice viscosity.
106    
107     @return: relative tolerance
108     @rtype rtol: positive C{float}
109     """
110     return self.__rtol
111     #===========================================================================
112     def setDruckerPragerLaw(self,tau_Y=None,friction=0):
113     """
114     Sets the parameters for the Drucker-Prager model.
115    
116     @param tau_Y: yield stress
117     @param friction: friction coefficient
118     """
119     self.__tau_Y=tau_Y
120     self.__friction=friction
121     def getFriction(self):
122     """
123     returns the friction coefficient
124    
125     @return: friction coefficient
126     """
127     return self.__friction
128     def getTauY(self):
129     """
130     returns the yield stress
131    
132     @return: the yield stress
133     """
134     return self.__tau_Y
135     #===========================================================================
136     def getElasticShearModulus(self):
137     """
138     returns the elastic shear modulus.
139    
140     @return: elastic shear modulus
141     """
142     return self.__mu
143     def setElasticShearModulus(self,mu=None):
144     """
145     Sets the elastic shear modulus.
146    
147     @param mu: elastic shear modulus
148     """
149     self.__mu=mu
150     #===========================================================================
151     def getPower(self, id=None):
152     """
153     returns the power in the power law
154    
155     @param id: if present, the power for material C{id} is returned.
156     @type id: C{int}
157     @return: the list of the powers for all matrials is returned. If C{id} is present only the power for material C{id} is returned.
158     """
159     if id == None:
160     return self.__power
161     else:
162     if self.validMaterialId(id):
163     return self.__power[id]
164     else:
165     raise ValueError,"Illegal material id %s."%id
166     def getEtaN(self, id=None):
167     """
168     returns the viscosity
169    
170     @param id: if present, the viscosity for material C{id} is returned.
171     @type id: C{int}
172     @return: the list of the viscosities for all matrials is returned. If C{id} is present only the viscosity for material C{id} is returned.
173     """
174     if id == None:
175     return self.__eta_N
176     else:
177     if self.validMaterialId(id):
178     return self.__eta_N[id]
179     else:
180     raise ValueError,"Illegal material id %s."%id
181     def getTauT(self, id=None):
182     """
183     returns the transition stress
184    
185     @param id: if present, the transition stress for material C{id} is returned.
186     @type id: C{int}
187     @return: the list of the transition stresses for all matrials is returned. If C{id} is present only the transition stress for material C{id} is returned.
188     """
189     if id == None:
190     return self.__tau_t
191     else:
192     if self.validMaterialId(id):
193     return self.__tau_t[id]
194     else:
195     raise ValueError,"Illegal material id %s."%id
196    
197     def setPowerLaw(self,eta_N, id=0, tau_t=1, power=1):
198     """
199     Sets the power-law parameters for material id
200    
201     @param id: material id
202     @type id: C{int}
203     @param eta_N: viscosity for tau=tau_t
204     @param tau_t: transition stress
205     @param power: power law coefficient
206     """
207     if self.validMaterialId(id):
208     self.__eta_N[id]=eta_N
209     self.__power[id]=power
210     self.__tau_t[id]=tau_t
211     else:
212     raise ValueError,"Illegal material id %s."%id
213    
214     def setPowerLaws(self,eta_N, tau_t, power):
215     """
216     Sets the parameters of the power-law for all materials.
217    
218     @param eta_N: list of viscosities for tau=tau_t
219     @param tau_t: list of transition stresses
220     @param power: list of power law coefficient
221     """
222     if len(eta_N)!=self.__numMaterials or len(tau_t)!=self.__numMaterials or len(power)!=self.__numMaterials:
223     raise ValueError,"%s materials are expected."%self.__numMaterials
224     for i in xrange(self.__numMaterials):
225     self.setPowerLaw(id=i, eta_N=eta_N[i],tau_t=tau_t[i],power=power[i])
226    
227     #===========================================================================
228     def getEtaEff(self,gamma_dot, eta0=None, pressure=0,dt=None, iter_max=10):
229     """
230     returns the effective viscosity eta_eff such that
231    
232     M{tau=eta_eff * gamma_dot}
233    
234     by solving a non-linear problem for tau.
235    
236     @param gamma_dot: equivalent strain gamma_dot
237     @param eta0: initial guess for the effective viscosity (e.g from a previous time step). If not present, an initial guess is calculated.
238     @param pressure: pressure used to calculate yield condition
239     @param dt: time step size. only needed if elastic component is considered.
240     @type dt: positive C{float} if present
241     @param iter_max: maximum number of iteration steps.
242     @type iter_max: C{int}
243     @return: effective viscosity.
244     """
245     numMaterial=self.getNumMaterials()
246     s=[1.-1./p for p in self.getPower() ]
247     eta_N=self.getEtaN()
248     tau_t=self.getTauT()
249     mu=self.getElasticShearModulus()
250     fric=self.getFriction()
251     tau_Y=self.getTauY()
252     if eta0==None:
253     theta=0.
254     for i in xrange(numMaterial):
255     inv_eta_i=0**s[i]/eta_N[i]
256     theta=theta+inv_eta_i
257     if util.inf(theta)<=0:
258     raise ValueError,"unable to set positive initial guess for eta_eff. Most likely no power law with power 1 set."
259     eta_eff=1./theta
260     else:
261     if util.inf(eta0)<=0:
262     raise ValueError,"initial guess for eta_eff is not positive."
263     eta_eff=eta0
264    
265     if mu !=None and dt == None:
266     raise ValueError,"Time stepsize dt must be given."
267     if dt !=None:
268     if dt<=0: raise ValueEror,"time step size must be positive."
269     if tau_Y==None:
270     eta_max=None
271     else:
272     if util.inf(pressure)*fric<0:
273     raise ValueError,"pressure needs to be non-negative."
274     eta_max=(tau_Y+fric*pressure)/(gamma_dot+util.whereZero(gamma_dot)/util.DBLE_MAX)
275     rtol=self.getEtaTolerance()
276     iter =0
277     converged=False
278     tau=eta_eff*gamma_dot
279     if self.__verbose: print "Start calculation of eta_eff (tolerance = %s)\ninitial max eta_eff = %s, tau = %s."%(rtol,util.Lsup(eta_eff),util.Lsup(tau))
280     while not converged:
281     if iter>max(iter_max,1):
282     raise RunTimeError,"tolerance not reached after %s steps."%max(iter_max,1)
283     #===========================================
284     theta=0. # =1/eta
285     omega=0. # = tau*theta'= eta'*tau/eta**2
286     if mu !=None: theta=1./(dt*mu)
287     for i in xrange(numMaterial):
288     inv_eta_i=(tau/tau_t[i])**s[i]/eta_N[i]
289     theta=theta+inv_eta_i
290     omega=omega+s[i]*inv_eta_i
291     #===========================================
292     eta_eff, eta_eff_old=util.clip(eta_eff*(theta+omega)/(eta_eff*theta**2+omega),maxval=eta_max), eta_eff
293     tau=eta_eff*gamma_dot
294     d=util.Lsup(eta_eff-eta_eff_old)
295     l=util.Lsup(eta_eff)
296     iter+=1
297     if self.__verbose: print "step %s: correction = %s, max eta_eff = %s, max tau= %s"%(iter, d, l,util.Lsup(tau))
298     converged= d<= rtol* l
299     return eta_eff
300    
301     #====================================================================================================================================
302    
303 gross 2100 class IncompressibleIsotropicKelvinFlow(Defect):
304 gross 1639 """
305 caltinay 2169 This class implements the rheology of an isotropic Kelvin material.
306 gross 2100
307 caltinay 2169 Typical usage::
308 gross 1639
309 caltinay 2169 sp = IncompressibleIsotropicKelvinFlow(domain, stress=0, v=0)
310     sp.setTolerance()
311     sp.initialize(...)
312     v,p = sp.solve(v0, p0)
313 gross 1639
314 caltinay 2169 @note: This model has been used in the self-consistent plate-mantle model
315     proposed in U{Hans-Bernd Muhlhaus<emailto:h.muhlhaus@uq.edu.au>}
316     and U{Klaus Regenauer-Lieb<mailto:klaus.regenauer-lieb@csiro.au>}:
317     I{Towards a self-consistent plate mantle model that includes elasticity: simple benchmarks and application to basic modes of convection},
318     see U{doi: 10.1111/j.1365-246X.2005.02742.x<http://www3.interscience.wiley.com/journal/118661486/abstract?CRETRY=1&SRETRY=0>}.
319    
320 gross 1639 """
321 gross 2100 def __init__(self, domain, stress=0, v=0, p=0, t=0, numMaterials=1, useJaumannStress=True, **kwargs):
322 gross 1639 """
323 caltinay 2169 Initializes the model.
324 gross 1639
325     @param domain: problem domain
326     @type domain: L{domain}
327 gross 2100 @param stress: initial deviatoric stress
328 gross 1639 @param v: initial velocity field
329 gross 2100 @param p: initial pressure
330 gross 1639 @param t: initial time
331 caltinay 2169 @param useJaumannStress: C{True} if Jaumann stress is used
332     (not supported yet)
333 gross 1639 """
334 gross 2100 super(IncompressibleIsotropicKelvinFlow, self).__init__(**kwargs)
335 gross 1639 self.__domain=domain
336     self.__t=t
337     self.__vol=util.integrate(1.,Function(self.__domain))
338     self.__useJaumannStress=useJaumannStress
339 gross 2100 self.__v_boundary=Vector(0,Solution(self.__domain))
340 gross 1639 #=======================
341     # state variables:
342     #
343     if isinstance(stress,Data):
344     self.__stress=util.interpolate(stress,Function(domain))
345     else:
346     self.__stress=Data(stress,(domain.getDim(),domain.getDim()),Function(domain))
347     self.__stress-=util.trace(self.__stress)*(util.kronecker(domain)/domain.getDim())
348     if isinstance(v,Data):
349     self.__v=util.interpolate(v,Solution(domain))
350     else:
351 caltinay 2169 self.__v=Data(v,(domain.getDim(),),Solution(domain))
352 gross 1639 if isinstance(p,Data):
353     self.__p=util.interpolate(p,ReducedSolution(domain))
354     else:
355     self.__p=Data(p,(),ReducedSolution(domain))
356     #=======================
357 caltinay 2169 # PDE related stuff
358 gross 2100 self.__pde_v=LinearPDE(domain,numEquations=self.getDomain().getDim(),numSolutions=self.getDomain().getDim())
359     self.__pde_v.setSymmetryOn()
360     self.__pde_v.setSolverMethod(preconditioner=LinearPDE.RILU)
361 caltinay 2169
362 gross 2100 self.__pde_p=LinearPDE(domain)
363     self.__pde_p.setReducedOrderOn()
364     self.__pde_p.setSymmetryOn()
365 gross 1639
366 gross 2100 self.setTolerance()
367     self.setSmall()
368 caltinay 2169
369 gross 1639 def useJaumannStress(self):
370     """
371 caltinay 2169 Returns C{True} if Jaumann stress is included.
372     """
373 gross 1639 return self.__useJaumannStress
374 gross 2100
375 gross 1639 def setSmall(self,small=util.sqrt(util.EPSILON)):
376     """
377 caltinay 2169 Sets a small value to be used.
378    
379 gross 1639 @param small: positive small value
380     """
381     self.__small=small
382 gross 2100
383 gross 1639 def getSmall(self):
384     """
385 caltinay 2169 Returns small value.
386 gross 1639 @rtype: positive float
387     """
388     return self.__small
389 gross 2100
390     def setTolerance(self,tol=1.e-4):
391 gross 1639 """
392 caltinay 2169 Sets the tolerance.
393 gross 1639 """
394 gross 2100 self.__pde_v.setTolerance(tol**2)
395     self.__pde_p.setTolerance(tol**2)
396     self.__tol=tol
397    
398     def getTolerance(self):
399 gross 1639 """
400 caltinay 2169 Returns the set tolerance.
401 gross 2100 @rtype: positive float
402 gross 1639 """
403 gross 2100 return self.__tol
404    
405     def getDomain(self):
406 gross 1639 """
407 caltinay 2169 Returns the domain.
408 gross 1639 """
409 gross 2100 return self.__domain
410 caltinay 2169
411 gross 1639 def getTime(self):
412     """
413 caltinay 2169 Returns current time.
414 gross 1639 """
415     return self.__t
416 gross 2100
417     def setExternals(self, F=None, f=None, q=None, v_boundary=None):
418 caltinay 2169 """
419     Sets externals.
420 gross 1639
421 caltinay 2169 @param F: external force
422 gross 2100 @param f: surface force
423 gross 1639 @param q: location of constraints
424     """
425 gross 2100 if F != None: self.__pde_v.setValue(Y=F)
426     if f != None: self.__pde_v.setValue(y=f)
427     if q != None: self.__pde_v.setValue(q=q)
428     if v_boundary != None: self.__v_boundary=v_boundary
429 gross 1639
430 gross 2100 def bilinearform(self,arg0,arg1):
431     s0=util.deviatoric(util.symmetric(util.grad(arg0[0])))
432     s1=util.deviatoric(util.symmetric(util.grad(arg1[0])))
433     # s0=util.interpolate(arg0[0],Function(self.getDomain()))
434     # s1=util.interpolate(arg1[0],Function(self.getDomain()))
435     p0=util.interpolate(arg0[1],Function(self.getDomain()))
436     p1=util.interpolate(arg1[1],Function(self.getDomain()))
437     a=util.integrate(self.__p_weight**2*util.inner(s0,s1))+util.integrate(p0*p1)
438     return a
439    
440     def getEtaEff(self,strain, pressure):
441 caltinay 2169 if self.__mu==None:
442 gross 2100 eps=util.length(strain)*util.sqrt(2)
443 gross 1639 else:
444 gross 2100 eps=util.length(strain+self.__stress/((2*self.__dt)*self.__mu))*util.sqrt(2)
445     p=util.interpolate(pressure,eps.getFunctionSpace())
446     if self.__tau_Y!= None:
447     tmp=self.__tau_Y+self.__friction*p
448     m=util.wherePositive(eps)*util.wherePositive(tmp)
449     eta_max=m*tmp/(eps+(1-m)*util.EPSILON)+(1-m)*util.DBLE_MAX
450 gross 1639 else:
451 gross 2100 eta_max=util.DBLE_MAX
452     # initial guess:
453     tau=util.length(self.__stress)/util.sqrt(2)
454     # start the iteration:
455     cc=0
456     TOL=1e-7
457     dtau=util.DBLE_MAX
458     print "tau = ", tau, "eps =",eps
459     while cc<10 and dtau>TOL*util.Lsup(tau):
460     eta_eff2,eta_eff_dash=self.evalEtaEff(tau,return_dash=True)
461     eta_eff=util.clip(eta_eff2-eta_eff_dash*tau/(1-eta_eff_dash*eps),maxval=eta_max)
462     tau, tau2=eta_eff*eps, tau
463     dtau=util.Lsup(tau2-tau)
464     print "step ",cc,dtau, util.Lsup(tau)
465     cc+=1
466     return eta_eff
467 gross 1639
468 gross 2100 def getEtaCharacteristic(self):
469     a=0
470     for i in xrange(self.__numMaterials):
471     a=a+1./self.__eta_N[i]
472     return 1/a
473 caltinay 2169
474 gross 2100 def evalEtaEff(self, tau, return_dash=False):
475     a=Scalar(0,tau.getFunctionSpace()) # =1/eta
476     if return_dash: a_dash=Scalar(0,tau.getFunctionSpace()) # =(1/eta)'
477     s=util.Lsup(tau)
478     if s>0:
479     m=util.wherePositive(tau)
480     tau2=s*util.EPSILON*(1.-m)+m*tau
481     for i in xrange(self.__numMaterials):
482     eta_N=self.__eta_N[i]
483     tau_t=self.__tau_t[i]
484     if tau_t==None:
485     a+=1./eta_N
486     else:
487     power=1.-1./self.__power[i]
488     c=1./(tau_t**power*eta_N)
489     a+=c*tau2**power
490     if return_dash: a_dash+=power*c*tau2**(power-1.)
491     else:
492     for i in xrange(self.__numMaterials):
493     eta_N=self.__eta_N[i]
494     power=1.-1./self.__power[i]
495     a+=util.whereZero(power)/eta_N
496     if self.__mu!=None: a+=1./(self.__dt*self.__mu)
497     out=1/a
498     if return_dash:
499     return out,-out**2*a_dash
500     else:
501     return out
502 caltinay 2169
503 gross 2100 def eval(self,arg):
504     v=arg[0]
505     p=arg[1]
506     D=self.getDeviatoricStrain(v)
507     eta_eff=self.getEtaEff(D,p)
508     print "eta_eff=",eta_eff
509     # solve for dv
510     self.__pde_v.setValue(A=Data()) # save memory!
511     k3=util.kronecker(Function(self.getDomain()))
512     k3Xk3=util.outer(k3,k3)
513     self.__pde_v.setValue(A=eta_eff*(util.swap_axes(k3Xk3,0,3)+util.swap_axes(k3Xk3,1,3)),X=-eta_eff*D+p*util.kronecker(self.getDomain()))
514     dv=self.__pde_v.getSolution(verbose=self.__verbose)
515     print "resistep dv =",dv
516     # solve for dp
517     v2=v+dv
518     self.__pde_p.setValue(D=1/eta_eff,Y=util.div(v2))
519     dp=self.__pde_p.getSolution(verbose=self.__verbose)
520     print "resistep dp =",dp
521     return ArithmeticTuple(dv,dp)
522    
523 caltinay 2169 def update(self, dt, iter_max=100, inner_iter_max=20, verbose=False):
524 gross 1639 """
525 caltinay 2169 Updates stress, velocity and pressure for time increment dt.
526 gross 1639
527     @param dt: time increment
528 caltinay 2169 @param inner_iter_max: maximum number of iteration steps in the
529     incompressible solver
530     @param verbose: prints some infos in the incompressible solver
531 gross 1639 """
532 gross 2100 self.__verbose=verbose
533     self.__dt=dt
534     tol=self.getTolerance()
535     # set the initial velocity:
536     m=util.wherePositive(self.__pde_v.getCoefficient("q"))
537     v_new=self.__v*(1-m)+self.__v_boundary*m
538     # and off we go:
539     x=ArithmeticTuple(v_new, self.__p)
540     # self.__p_weight=util.interpolate(1./self.getEtaCharacteristic(),Function(self.__domain))**2
541     self.__p_weight=self.getEtaCharacteristic()
542     # self.__p_weight=util.interpolate(1./self.getEtaCharacteristic()**2,self.__p.getFunctionSpace())
543     atol=self.norm(x)*self.__tol
544     x_new=NewtonGMRES(self, x, iter_max=iter_max,sub_iter_max=inner_iter_max, atol=atol,rtol=0., verbose=verbose)
545     self.__v=x_new[0]
546     self.__p=x_new[1]
547     1/0
548     # self.__stress=self.getUpdatedStress(...)
549     self.__t+=dt
550     return self.__v, self.__p
551 gross 1639
552 caltinay 2169 #========================================================================
553 gross 1639
554 gross 2100 def getNewDeviatoricStress(self,D,eta_eff=None):
555     if eta_eff==None: eta_eff=self.evalEtaEff(self.__stress,D,self.__p)
556     s=(2*eta_eff)*D
557     if self.__mu!=None: s+=eta_eff/(self.__dt*self.__mu)*self.__last_stress
558     return s
559 gross 1639
560 gross 2100 def getDeviatoricStress(self):
561     """
562 caltinay 2169 Returns current stress.
563 gross 2100 """
564     return self.__stress
565 caltinay 2169
566 gross 2100 def getDeviatoricStrain(self,velocity=None):
567     """
568 caltinay 2169 Returns strain.
569 gross 2100 """
570     if velocity==None: velocity=self.getVelocity()
571     return util.deviatoric(util.symmetric(util.grad(velocity)))
572 gross 1659
573 gross 2100 def getPressure(self):
574     """
575 caltinay 2169 Returns current pressure.
576 gross 2100 """
577     return self.__p
578 caltinay 2169
579 gross 2100 def getVelocity(self):
580     """
581 caltinay 2169 Returns current velocity.
582 gross 2100 """
583     return self.__v
584 gross 1639
585 gross 2100 def getTau(self,stress=None):
586 gross 1639 """
587 caltinay 2169 Returns current second stress deviatoric invariant.
588 gross 1639 """
589 gross 2100 if stress==None: stress=self.getDeviatoricStress()
590     return util.sqrt(0.5)*util.length(stress)
591 gross 1639
592 gross 2100 def getGammaDot(self,strain=None):
593     """
594 caltinay 2169 Returns current second stress deviatoric invariant.
595 gross 2100 """
596     if strain==None: strain=self.getDeviatoricStrain()
597     return util.sqrt(2)*util.length(strain)
598 gross 1639

  ViewVC Help
Powered by ViewVC 1.1.26