1 # -*- coding: utf-8 -*-
3 # Copyright (C) 2008-2024 EDF R&D
5 # This library is free software; you can redistribute it and/or
6 # modify it under the terms of the GNU Lesser General Public
7 # License as published by the Free Software Foundation; either
8 # version 2.1 of the License.
10 # This library is distributed in the hope that it will be useful,
11 # but WITHOUT ANY WARRANTY; without even the implied warranty of
12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
13 # Lesser General Public License for more details.
15 # You should have received a copy of the GNU Lesser General Public
16 # License along with this library; if not, write to the Free Software
17 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
19 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
21 # Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
23 #===============================================================================
25 # donnees numeriques du modele
32 p = 10 # nb de mesures
34 class EqChaleur(object):
36 def __init__(self,alpha):
38 self.alpha = float(alpha)
41 self.A=2*idx2*numpy.identity(dim)
42 for i in range(1,dim-2):
47 return numpy.linalg.solve(self.alpha*self.A,S)
51 S = numpy.zeros((dim,))
56 S[i] = pi2*numpy.sin(i*c)
59 def extract(self,champs,nbobs):
66 X = numpy.linspace(xdeb,xfin,dim)
69 return self.extract(T,nbobs)
73 X = numpy.linspace(xdeb,xfin,dim)
77 from matplotlib import pylab
81 class ExempleOperateurChaleur:
83 Modelisation d'operateur non lineaire pour la diffusion de la chaleur dans
84 un mur d'epaisseur e, avec une production de chaleur r, pour un
85 environnement a temperature Ts. Les equations sont :
87 T = Ts - r/2k x^2 + r e /2k x
89 def __init__(self, e = 10., Ts = 100., nbpt = 1000, nbobs = 10):
93 self.dx = self.e/(self.nbpt-1)
94 self.nbobs = min(int(nbobs),self.nbpt)
96 def Direct(self, XX ):
97 if type(XX) is type(numpy.matrix([])): alpha = XX.A1
98 elif type(XX) is type(numpy.array([])): alpha = numpy.matrix(XX).A1
102 HX = eq.obs(self.nbobs)
104 return numpy.array( HX )
106 def DirectExp(self, XX ): # Version en Log/Exp
107 if type(XX) is type(numpy.matrix([])): alpha = XX.A1
108 elif type(XX) is type(numpy.array([])): alpha = numpy.matrix(XX).A1
110 alpha = numpy.exp(alpha)
113 HX = eq.obs(self.nbobs)
115 return numpy.array( HX )
118 def Analytique(self, XX ):
121 for i in range(self.nbobs+1) :
122 ex = i * self.e/(self.nbobs)
123 eT = self.Ts - amort * ex**2 + self.e*amort*ex
124 listresult.append(eT)
125 return numpy.array(listresult)
127 def Tangent(self, paire ):
131 def Adjoint(self, paire ):
135 #===============================================================================
139 AlgorithmParameters = {
140 "EpsilonMinimumExponent" : -10,
141 # "PlotAndSave" : True,
143 "InitialDirection":CheckingPoint,
144 "AmplitudeOfInitialDirection":0.5,
148 OP = ExempleOperateurChaleur(nbpt = taille)
150 DirectOperator = OP.Direct
152 # JPA : probleme : il faut que le nom d'origine, qui est une méthode de classe,
153 # existe explicitement pour permettre le parallelisme !
154 Direct = DirectOperator