src/daComposant/daNumerics/Models/TwoDimensionalInverseDistanceCS2010.py

   1 # -*- coding: utf-8 -*-
   2 #
   3 # Copyright (C) 2008-2024 EDF R&D
   4 #
   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.
   9 #
  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.
  14 #
  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
  18 #
  19 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  20 #
  21 # Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
  22
  23 import sys, unittest, numpy
  24
  25 # ==============================================================================
  26 class TwoDimensionalInverseDistanceCS2010:
  27     """
  28     Two-dimensional inverse distance function of parameter µ=(µ1,µ2):
  29
  30         s(x,y,µ) = 1 / sqrt( (x-µ1)² + (y-µ2)² + 0.1² )
  31
  32         with (x,y) ∈ [0.1,0.9]² and µ=(µ1,µ2) ∈ [-1,-0.01]².
  33
  34     This is the non-linear parametric function (3.38) of the reference:
  35         Chaturantabut, S., Sorensen, D. C.,
  36         Nonlinear Model Reduction via Discrete Empirical Interpolation,
  37         SIAM Journal on Scientific Computing, 32(5), pp. 2737-2764 (2010).
  38     """
  39     def __init__(self, nx: int = 20, ny: int = 20):
  40         "Définition du maillage spatial"
  41         self.nx  = max(1, nx)
  42         self.ny  = max(1, ny)
  43         self.x   = numpy.linspace(0.1, 0.9, self.nx, dtype=float)
  44         self.y   = numpy.linspace(0.1, 0.9, self.ny, dtype=float)
  45
  46     def FieldG(self, mu ):
  47         "Fonction simulation pour un paramètre donné"
  48         mu1, mu2 = numpy.ravel( mu )
  49         #
  50         x, y = numpy.meshgrid( self.x, self.y )
  51         sxymu = 1. / numpy.sqrt( (x - mu1)**2 + (y - mu2)**2 + 0.1**2 )
  52         #
  53         return sxymu
  54
  55     def get_x(self):
  56         "Renvoie le maillage spatial"
  57         return self.x, self.y
  58
  59     def get_sample_of_mu(self, ns1: int = 20, ns2: int = 20):
  60         "Renvoie l'échantillonnage paramétrique régulier"
  61         smu1 = numpy.linspace(-1, -0.01, ns1, dtype=float)
  62         smu2 = numpy.linspace(-1, -0.01, ns2, dtype=float)
  63         smu = numpy.array([(mu1, mu2) for mu1 in smu1 for mu2 in smu2])
  64         return smu
  65
  66     def get_random_sample_of_mu(self, ns1: int = 1, ns2: int = 1):
  67         "Renvoie l'échantillonnage paramétrique aléatoire"
  68         smu = []
  69         for i in range(ns1 * ns2):
  70             smu1 = numpy.random.uniform(-1, -0.01)
  71             smu2 = numpy.random.uniform(-1, -0.01)
  72             smu.append((smu1, smu2))
  73         smu = numpy.array(smu)
  74         return smu
  75
  76     def get_bounds_on_space(self):
  77         "Renvoie les bornes sur le maillage spatial"
  78         return [[min(self.x), max(self.x)], [min(self.y), max(self.y)]]
  79
  80     def get_bounds_on_parameter(self):
  81         "Renvoie les bornes sur le maillage paramétrique"
  82         return [[-1, -0.01]] * 2
  83
  84     OneRealisation = FieldG
  85
  86 # ==============================================================================
  87 class LocalTest(unittest.TestCase):
  88     @classmethod
  89     def setUpClass(cls):
  90         print('\nAUTODIAGNOSTIC\n==============\n')
  91         print("    " + TwoDimensionalInverseDistanceCS2010().__doc__.strip())
  92
  93     def test001(self):
  94         numpy.random.seed(123456789)
  95         Equation = TwoDimensionalInverseDistanceCS2010()
  96         for mu in Equation.get_sample_of_mu(5, 5):
  97             solution = Equation.OneRealisation( mu )
  98             # Nappe maximale au coin (0,0)
  99             self.assertTrue(numpy.max(solution.flat) <= solution[0, 0])
 100             # Nappe minimale au coin [-1,-1]
 101             self.assertTrue(numpy.min(solution.flat) >= solution[-1, -1])
 102
 103     def tearDown(cls):
 104         print("\n    Tests OK\n")
 105
 106 # ==============================================================================
 107 if __name__ == "__main__":
 108     sys.stderr = sys.stdout
 109     unittest.main(verbosity=0)