.. include:: snippets/Header2Algo01.rst
This algorithm provides a simple way of analyzing the quality of the empirical
-state interpolation obtained using a reduced basis, using measurements at
-specific points.
+interpolation obtained by a reduced basis for complete states, using
+measurements at precise points.
The results displayed by default are simple statistics related to the
-normalized errors of interpolation with a reduced basis. The test uses a
-reduced base and a set of optimal measurement positions, and uses
-pseudo-measurements from each complete state (“*snapshot*”) included in the
+normalized errors between the interpolation with a reduced basis and the
+complete states we seek to represent on the reduced basis. The test requires a
+reduced basis and a set of optimal measurement positions, and uses
+pseudo-measurements from each complete state ("*snapshot*") included in the
given test set.
Please note: to be consistent, this test must use the same mathematical norm as
.. literalinclude:: scripts/simple_MeasurementsOptimalPositioningTask2.res
:language: none
+.. --------- ..
+.. include:: scripts/simple_MeasurementsOptimalPositioningTask3.rst
+
+.. literalinclude:: scripts/simple_MeasurementsOptimalPositioningTask3.py
+
+.. include:: snippets/Header2Algo10.rst
+
+.. literalinclude:: scripts/simple_MeasurementsOptimalPositioningTask3.res
+ :language: none
+
+.. include:: snippets/Header2Algo11.rst
+
+.. _simple_MeasurementsOptimalPositioningTask31:
+.. image:: scripts/simple_MeasurementsOptimalPositioningTask31.png
+ :align: center
+ :width: 60%
+
+.. _simple_MeasurementsOptimalPositioningTask32:
+.. image:: scripts/simple_MeasurementsOptimalPositioningTask32.png
+ :align: center
+ :width: 60%
+
+.. _simple_MeasurementsOptimalPositioningTask33:
+.. image:: scripts/simple_MeasurementsOptimalPositioningTask33.png
+ :align: center
+ :width: 60%
+
.. ------------------------------------ ..
.. include:: snippets/Header2Algo06.rst
- data text file (TXT, CSV, TSV, DAT), with variable pointer by name in column or row, indicated by the keyword "*DataFile*" with the condition ``Vector=True``
- binary data file (NPY, NPZ), with variable pointer by name, indicated by the keyword "*DataFile*" with the condition ``Vector=True``
-- Examples of statements in TUI interface:
+Examples of statements in TUI interface:
+::
- - ``case.setObservation( Vector = [1, 2, 3] )``
- - ``case.setObservation( Vector = numpy.array([1, 2, 3]) )``
- - ``case.setObservation( Vector = '1 2 3' )``
- - ``case.setObservation( Vector=True, Script = 'script.py' )```
- - ``case.setObservation( Vector=True, DataFile = 'data.csv' )```
- - ``case.setObservation( Vector=True, DataFile = 'data.npy' )```
+ case.setObservation( Vector = [1, 2, 3] )
+ case.setObservation( Vector = numpy.array([1, 2, 3]) )
+ case.setObservation( Vector = '1 2 3' )
+ case.setObservation( Vector=True, Script = 'script.py' )
+ case.setObservation( Vector=True, DataFile = 'data.csv' )
+ case.setObservation( Vector=True, DataFile = 'data.npy' )
Use note: in a given study, only the last record (whether a single vector or a
series of vectors) can be used, as only one observation concept exists per ADAO
- data text file (TXT, CSV, TSV), with variable pointer by name in column or row, indicated by the keyword "*DataFile*" with the condition ``VectorSerie=True``
- binary data file (NPY, NPZ), with variable pointer by name, indicated by the keyword "*DataFile*" with the condition ``VectorSerie=True``
-- Examples of statements in TUI interface:
-
- - ``case.setObservation( VectorSerie = [[1,2,3], [1,2,3]] )``
- - ``case.setObservation( VectorSerie = [numpy.array([1,2,3]), numpy.array([1,2,3])] )``
- - ``case.setObservation( VectorSerie = ['1 2 3', '1 2 3'] )``
- - ``case.setObservation( VectorSerie = '[[1,2,3], [1,2,3]]' )``
- - ``case.setObservation( VectorSerie = '1 2 3 ; 1 2 3' )``
- - ``case.setObservation( VectorSerie=True, Script = 'script.py' )```
- - ``case.setObservation( VectorSerie=True, DataFile = 'data.csv' )```
- - ``case.setObservation( VectorSerie=True, DataFile = 'data.npy' )```
+Examples of statements in TUI interface:
+::
+
+ case.setObservation( VectorSerie = [[1,2,3], [1,2,3]] )
+ case.setObservation( VectorSerie = [numpy.array([1,2,3]), numpy.array([1,2,3])] )
+ case.setObservation( VectorSerie = ['1 2 3', '1 2 3'] )
+ case.setObservation( VectorSerie = '[[1,2,3], [1,2,3]]' )
+ case.setObservation( VectorSerie = '1 2 3 ; 1 2 3' )
+ case.setObservation( VectorSerie=True, Script = 'script.py' )
+ case.setObservation( VectorSerie=True, DataFile = 'data.csv' )
+ case.setObservation( VectorSerie=True, DataFile = 'data.npy' )
Use note: in a given study, only the last record (whether a single vector or a
series of vectors) can be used, as only one observation concept exists per ADAO
--- /dev/null
+# -*- coding: utf-8 -*-
+#
+import numpy as np
+import matplotlib.pyplot as plt
+np.random.seed(123456789)
+#
+dimension = 20
+#
+print("Defining a set of artificial physical fields")
+from Models.TwoDimensionalInverseDistanceCS2010 \
+ import TwoDimensionalInverseDistanceCS2010 as Equation
+Eq = Equation(dimension, dimension)
+print()
+#
+print("Search for optimal measurement positions")
+from adao import adaoBuilder
+case = adaoBuilder.New()
+case.setAlgorithmParameters(
+ Algorithm = 'MeasurementsOptimalPositioningTask',
+ Parameters = {
+ "Variant":"DEIM",
+ "SampleAsnUplet":Eq.get_sample_of_mu(15, 15),
+ "MaximumNumberOfLocations":50,
+ "ErrorNorm":"Linf",
+ "ErrorNormTolerance":0.,
+ "StoreSupplementaryCalculations":[
+ "SingularValues",
+ "OptimalPoints",
+ "EnsembleOfSimulations",
+ ],
+ }
+ )
+case.setBackground(Vector = [1,1] )
+case.setObservationOperator(OneFunction = Eq.OneRealisation)
+case.execute()
+#
+print()
+print("Graphical display of the results")
+#
+sp = case.get("EnsembleOfSimulations")[-1]
+x1, x2 = Eq.get_x()
+x1, x2 = np.meshgrid(x1, x2)
+name = "Representation of a few snapshots of G out of a total of %i"%sp.shape[1]
+print(" -", name)
+fig, ax = plt.subplots(subplot_kw={"projection": "3d"}, figsize=(6, 6))
+fig.suptitle(name)
+ax.set_xlabel("Position x1", fontweight='bold', color='red')
+ax.set_ylabel("Position x2", fontweight='bold', color='red')
+ax.set_zlabel("Magnitude of G", fontweight='bold', color='red')
+for i in range(sp.shape[1]):
+ if i % 44 != 0: continue
+ ax.plot_surface(x1, x2, sp[:,i].reshape((dimension, dimension)), cmap='coolwarm')
+fig.savefig("simple_MeasurementsOptimalPositioningTask31.png")
+plt.close()
+#
+sv = case.get("SingularValues")[-1]
+name = "Singular values of the set of G simulations"
+print(" -", name)
+fig, ax = plt.subplots(figsize=(6, 6))
+fig.suptitle(name)
+ax.set_xlabel("Index of singular values, numbered from 1")
+ax.set_ylabel("Magnitude of singular values")
+ax.set_xlim(1, len(sv))
+ax.set_yscale("log")
+ax.grid(True)
+ax.plot(range(1, 1 + len(sv)), sv)
+fig.savefig("simple_MeasurementsOptimalPositioningTask32.png")
+plt.tight_layout()
+plt.close()
+#
+nbmax = 15
+op = case.get("OptimalPoints")[-1]
+posx1 = [x1.reshape((-1,))[ip] for ip in op[:nbmax]]
+posx2 = [x2.reshape((-1,))[ip] for ip in op[:nbmax]]
+name = "Set of %i first optimal points of measurement"%nbmax
+Omega = Eq.get_bounds_on_space()
+print(" -", name)
+fig, ax = plt.subplots(figsize=(6, 6))
+fig.suptitle(name)
+ax.set_xlabel("Position x1", fontweight='bold', color='red')
+ax.set_ylabel("Position x2", fontweight='bold', color='red')
+ax.set_xlim(Omega[0][0] - 0.05, Omega[0][1] + 0.05)
+ax.set_ylim(Omega[1][0] - 0.05, Omega[1][1] + 0.05)
+ax.grid(True, which='both', linestyle=(0, (1, 5)), linewidth=0.5)
+ax.plot(posx1, posx2, markersize=6, marker="o", linestyle='')
+for i in range(len(posx1)):
+ ax.text(posx1[i] + 0.005, posx2[i] + 0.01, str(i + 1), fontweight='bold')
+fig.savefig("simple_MeasurementsOptimalPositioningTask33.png")
+plt.tight_layout()
+plt.close()
--- /dev/null
+Defining a set of artificial physical fields
+
+Search for optimal measurement positions
+
+Graphical display of the results
+ - Representation of a few snapshots of G out of a total of 225
+ - Singular values of the set of G simulations
+ - Set of 15 first optimal points of measurement
--- /dev/null
+.. index:: single: MeasurementsOptimalPositioningTask (example)
+
+Third example
+.............
+
+This more complete example describes the optimal positioning of measurements
+associated with a reduced DEIM-type decomposition on a classical parametric
+non-linear function, as proposed in reference [Chaturantabut10]_. This
+particular function has the notable advantage of being dependent only on a
+position :math:`x` in 2D and a parameter :math:`\mu` in 2D too. It therefore
+enables a pedagogical illustration of optimal measurement points.
+
+This function depends on the position :math:`x=(x_1,x_2)\in\Omega=[0.1,0.9]^2`
+in the 2D plane, and on the parameter
+:math:`\mu=(\mu_1,\mu_2)\in\mathcal{D}=[-1,-0.01]^2` of dimension 2 :
+
+.. math:: G(x;\mu) = \frac{1}{\sqrt{(x_1 - \mu_1)^2 + (x_2 - \mu_2)^2 + 0.1^2}}
+
+The function is represented on a regular :math:`\Omega_G` spatial grid of size
+20x20 points. It is available in ADAO built-in test models under the name
+`TwoDimensionalInverseDistanceCS2010`. So here we first build a set of
+simulations of :math:`G`, then apply the DEIM-type decomposition algorithm to
+it, and derive some simple illustrations.
+
+It can be seen that the singular values decrease steadily down to numerical
+noise, indicating that around a hundred basis elements are needed to fully
+represent the information contained in the set of :math:`G` simulations.
+Furthermore, the optimal measurement points in the :math:`\Omega_G` domain are
+inhomogeneously distributed, favoring the spatial zone near the
+:math:`(0.1,0.1)` corner in which the :math:`G` function varies more.
.. include:: snippets/Header2Algo01.rst
Cet algorithme permet d'analyser de manière simple la qualité de
-l'interpolation empirique des états obtenue par une base réduite, en utilisant
-des mesures en des points précis.
+l'interpolation empirique obtenue par une base réduite pour des états complets,
+en utilisant des mesures en des points précis.
Les résultats affichés par défaut sont des statistiques simples liées aux
-erreurs normalisées de l'interpolation avec une base réduite. Le test utilise
+erreurs normalisées entre l'interpolation avec une base réduite et les états
+complets que l'on cherche à représenter sur la base réduite. Le test requiert
une base réduite et un ensemble de positions de mesures optimales, et utilise
des pseudo-mesures provenant de chaque état complet ("*snapshot*") inclus dans
-l'ensemble de test donné.
+l'ensemble donné de test.
Attention : pour être cohérent, ce test doit utiliser la même norme
mathématique que celle utilisée pour construire la base réduite. La norme étant
.. literalinclude:: scripts/simple_MeasurementsOptimalPositioningTask2.res
:language: none
+.. --------- ..
+.. include:: scripts/simple_MeasurementsOptimalPositioningTask3.rst
+
+.. literalinclude:: scripts/simple_MeasurementsOptimalPositioningTask3.py
+
+.. include:: snippets/Header2Algo10.rst
+
+.. literalinclude:: scripts/simple_MeasurementsOptimalPositioningTask3.res
+ :language: none
+
+.. include:: snippets/Header2Algo11.rst
+
+.. _simple_MeasurementsOptimalPositioningTask31:
+.. image:: scripts/simple_MeasurementsOptimalPositioningTask31.png
+ :align: center
+ :width: 60%
+
+.. _simple_MeasurementsOptimalPositioningTask32:
+.. image:: scripts/simple_MeasurementsOptimalPositioningTask32.png
+ :align: center
+ :width: 60%
+
+.. _simple_MeasurementsOptimalPositioningTask33:
+.. image:: scripts/simple_MeasurementsOptimalPositioningTask33.png
+ :align: center
+ :width: 60%
+
.. ------------------------------------ ..
.. include:: snippets/Header2Algo06.rst
- fichier texte de données (TXT, CSV, TSV, DAT), avec pointeur de variable par nom en colonne ou en ligne, indiqué par le mot-clé "*DataFile*" avec la condition ``Vector=True``
- fichier binaire de données (NPY, NPZ), avec pointeur de variable par nom, indiqué par le mot-clé "*DataFile*" avec la condition ``Vector=True``
-- Exemples de déclaration en interface TUI :
+Exemples de déclaration en interface TUI :
+::
- - ``case.setObservation( Vector = [1, 2, 3] )``
- - ``case.setObservation( Vector = numpy.array([1, 2, 3]) )``
- - ``case.setObservation( Vector = '1 2 3' )``
- - ``case.setObservation( Vector=True, Script = 'script.py' )```
- - ``case.setObservation( Vector=True, DataFile = 'data.csv' )```
- - ``case.setObservation( Vector=True, DataFile = 'data.npy' )```
+ case.setObservation( Vector = [1, 2, 3] )
+ case.setObservation( Vector = numpy.array([1, 2, 3]) )
+ case.setObservation( Vector = '1 2 3' )
+ case.setObservation( Vector=True, Script = 'script.py' )
+ case.setObservation( Vector=True, DataFile = 'data.csv' )
+ case.setObservation( Vector=True, DataFile = 'data.npy' )
Remarque d'utilisation : dans une étude donnée, seul le dernier enregistrement
(que ce soit un vecteur unique ou une série de vecteurs) est utilisable, car un
- fichier texte de données (TXT, CSV, TSV, DAT), avec pointeur de variable par nom en colonne ou en ligne, indiqué par le mot-clé "*DataFile*" avec la condition ``VectorSerie=True``
- fichier binaire de données (NPY, NPZ), avec pointeur de variable par nom, indiqué par le mot-clé "*DataFile*" avec la condition ``VectorSerie=True``
-- Exemples de déclaration en interface TUI :
-
- - ``case.setObservation( VectorSerie = [[1,2,3], [1,2,3]] )``
- - ``case.setObservation( VectorSerie = [numpy.array([1,2,3]), numpy.array([1,2,3])] )``
- - ``case.setObservation( VectorSerie = ['1 2 3', '1 2 3'] )``
- - ``case.setObservation( VectorSerie = '[[1,2,3], [1,2,3]]' )``
- - ``case.setObservation( VectorSerie = '1 2 3 ; 1 2 3' )``
- - ``case.setObservation( VectorSerie=True, Script = 'script.py' )```
- - ``case.setObservation( VectorSerie=True, DataFile = 'data.csv' )```
- - ``case.setObservation( VectorSerie=True, DataFile = 'data.npy' )```
+Exemples de déclaration en interface TUI :
+::
+
+ case.setObservation( VectorSerie = [[1,2,3], [1,2,3]] )
+ case.setObservation( VectorSerie = [numpy.array([1,2,3]), numpy.array([1,2,3])] )
+ case.setObservation( VectorSerie = ['1 2 3', '1 2 3'] )
+ case.setObservation( VectorSerie = '[[1,2,3], [1,2,3]]' )
+ case.setObservation( VectorSerie = '1 2 3 ; 1 2 3' )
+ case.setObservation( VectorSerie=True, Script = 'script.py' )
+ case.setObservation( VectorSerie=True, DataFile = 'data.csv' )
+ case.setObservation( VectorSerie=True, DataFile = 'data.npy' )
Remarque d'utilisation : dans une étude donnée, seul le dernier enregistrement
(que ce soit un vecteur unique ou une série de vecteurs) est utilisable, car un
--- /dev/null
+# -*- coding: utf-8 -*-
+#
+import numpy as np
+import matplotlib.pyplot as plt
+np.random.seed(123456789)
+#
+dimension = 20
+#
+print("Définition d'un ensemble artificiel de champs physiques")
+from Models.TwoDimensionalInverseDistanceCS2010 \
+ import TwoDimensionalInverseDistanceCS2010 as Equation
+Eq = Equation(dimension, dimension)
+print()
+#
+print("Recherche des positions optimales de mesure")
+from adao import adaoBuilder
+case = adaoBuilder.New()
+case.setAlgorithmParameters(
+ Algorithm = 'MeasurementsOptimalPositioningTask',
+ Parameters = {
+ "Variant":"DEIM",
+ "SampleAsnUplet":Eq.get_sample_of_mu(15, 15),
+ "MaximumNumberOfLocations":50,
+ "ErrorNorm":"Linf",
+ "ErrorNormTolerance":0.,
+ "StoreSupplementaryCalculations":[
+ "SingularValues",
+ "OptimalPoints",
+ "EnsembleOfSimulations",
+ ],
+ }
+ )
+case.setBackground(Vector = [1,1] )
+case.setObservationOperator(OneFunction = Eq.OneRealisation)
+case.execute()
+#
+print()
+print("Affichage graphique des résultats")
+#
+sp = case.get("EnsembleOfSimulations")[-1]
+x1, x2 = Eq.get_x()
+x1, x2 = np.meshgrid(x1, x2)
+name = "Représentation de quelques snapshots de G sur un total de %i"%sp.shape[1]
+print(" -", name)
+fig, ax = plt.subplots(subplot_kw={"projection": "3d"}, figsize=(6, 6))
+fig.suptitle(name)
+ax.set_xlabel("Position x1", fontweight='bold', color='red')
+ax.set_ylabel("Position x2", fontweight='bold', color='red')
+ax.set_zlabel("Amplitude de G", fontweight='bold', color='red')
+for i in range(sp.shape[1]):
+ if i % 44 != 0: continue
+ ax.plot_surface(x1, x2, sp[:,i].reshape((dimension, dimension)), cmap='coolwarm')
+fig.savefig("simple_MeasurementsOptimalPositioningTask31.png")
+plt.close()
+#
+sv = case.get("SingularValues")[-1]
+name = "Valeurs singulières de l'ensemble des simulations de G"
+print(" -", name)
+fig, ax = plt.subplots(figsize=(6, 6))
+fig.suptitle(name)
+ax.set_xlabel("Index des valeurs singulières, numérotées à partir de 1")
+ax.set_ylabel("Amplitude des valeurs singulières")
+ax.set_xlim(1, len(sv))
+ax.set_yscale("log")
+ax.grid(True)
+ax.plot(range(1, 1 + len(sv)), sv)
+fig.savefig("simple_MeasurementsOptimalPositioningTask32.png")
+plt.tight_layout()
+plt.close()
+#
+nbmax = 15
+op = case.get("OptimalPoints")[-1]
+posx1 = [x1.reshape((-1,))[ip] for ip in op[:nbmax]]
+posx2 = [x2.reshape((-1,))[ip] for ip in op[:nbmax]]
+name = "Ensemble des %i premiers points optimaux de mesure"%nbmax
+Omega = Eq.get_bounds_on_space()
+print(" -", name)
+fig, ax = plt.subplots(figsize=(6, 6))
+fig.suptitle(name)
+ax.set_xlabel("Position x1", fontweight='bold', color='red')
+ax.set_ylabel("Position x2", fontweight='bold', color='red')
+ax.set_xlim(Omega[0][0] - 0.05, Omega[0][1] + 0.05)
+ax.set_ylim(Omega[1][0] - 0.05, Omega[1][1] + 0.05)
+ax.grid(True, which='both', linestyle=(0, (1, 5)), linewidth=0.5)
+ax.plot(posx1, posx2, markersize=6, marker="o", linestyle='')
+for i in range(len(posx1)):
+ ax.text(posx1[i] + 0.005, posx2[i] + 0.01, str(i + 1), fontweight='bold')
+fig.savefig("simple_MeasurementsOptimalPositioningTask33.png")
+plt.tight_layout()
+plt.close()
--- /dev/null
+Définition d'un ensemble artificiel de champs physiques
+
+Recherche des positions optimales de mesure
+
+Affichage graphique des résultats
+ - Représentation de quelques snapshots de G sur un total de 225
+ - Valeurs singulières de l'ensemble des simulations de G
+ - Ensemble des 15 premiers points optimaux de mesure
--- /dev/null
+.. index:: single: MeasurementsOptimalPositioningTask (exemple)
+
+Troisième exemple
+.................
+
+Cet exemple plus complet décrit le positionnement optimal de mesures associé à
+une décomposition réduite de type DEIM sur une fonction classique paramétrique
+non-linéaire, telle que proposée dans la référence [Chaturantabut10]_. Cette
+fonction particulière a l'avantage notable de n'être dépendante que d'une
+position :math:`x` en 2D et d'un paramètre :math:`\mu` en 2D aussi. Elle permet
+donc une illustration pédagogique pour représenter les points optimaux de
+mesure.
+
+Cette fonction dépend de la position :math:`x=(x_1,x_2)\in\Omega=[0.1,0.9]^2`
+dans le plan 2D, et du paramètre
+:math:`\mu=(\mu_1,\mu_2)\in\mathcal{D}=[-1,-0.01]^2` de dimension 2 :
+
+.. math:: G(x;\mu) = \frac{1}{\sqrt{(x_1 - \mu_1)^2 + (x_2 - \mu_2)^2 + 0.1^2}}
+
+La fonction est représenté sur une grille spatiale régulière :math:`\Omega_G`
+de taille 20x20 points. Elle est disponible dans les modèles de tests intégrés
+pour ADAO sous le nom `TwoDimensionalInverseDistanceCS2010`. On construit donc
+ici tout d'abord un ensemble de simulations de :math:`G`, pour lui appliquer
+ensuite l'algorithme de décomposition de type DEIM, et en tirer des
+illustrations simples.
+
+On observe ainsi que les valeurs singulières décroissent régulièrement jusqu'au
+bruit numérique, indiquant qu'il faut environ une centaine d'éléments de base
+pour complètement représenter l'information contenue dans l'ensemble des
+simulations de :math:`G`. Par ailleurs, les points de mesure optimaux dans le
+domaine :math:`\Omega_G` sont répartis de manière inhomogène, privilégiant la
+zone spatiale proche du coin :math:`(0.1,0.1)` dans laquelle la fonction
+:math:`G` varie plus.
--- /dev/null
+# -*- coding: utf-8 -*-
+#
+# Copyright (C) 2008-2024 EDF R&D
+#
+# This library is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation; either
+# version 2.1 of the License.
+#
+# This library is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with this library; if not, write to the Free Software
+# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+#
+# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+# Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
+
+import sys, unittest, numpy
+
+# ==============================================================================
+class TwoDimensionalInverseDistanceCS2010:
+ """
+ Two-dimensional inverse distance function of parameter µ=(µ1,µ2):
+
+ s(x,y,µ) = 1 / sqrt( (x-µ1)² + (y-µ2)² + 0.1² )
+
+ avec (x,y) ∈ [0.1,0.9]² et µ=(µ1,µ2) ∈ [-1,-0.01]².
+
+ This is the non-linear parametric function (3.38) of the reference:
+ Chaturantabut, S., Sorensen, D. C.,
+ Nonlinear Model Reduction via Discrete Empirical Interpolation,
+ SIAM Journal on Scientific Computing, 32(5), pp. 2737-2764 (2010).
+ """
+ def __init__(self, nx: int = 20, ny: int = 20):
+ self.nx = max(1, nx)
+ self.ny = max(1, ny)
+ self.x = numpy.linspace(0.1, 0.9, self.nx, dtype=float)
+ self.y = numpy.linspace(0.1, 0.9, self.ny, dtype=float)
+
+ def FunctionH(self, XX ):
+ __mu1, __mu2 = numpy.ravel(XX)
+ #
+ __x, __y = numpy.meshgrid( self.x, self.y )
+ __sxymu = 1. / numpy.sqrt( (__x - __mu1)**2 + (__y - __mu2)**2 + 0.1**2 )
+ #
+ return __sxymu
+
+ def get_x(self):
+ return self.x, self.y
+
+ def get_sample_of_mu(self, ns1: int = 20, ns2: int = 20):
+ smu1 = numpy.linspace(-1, -0.01, ns1, dtype=float)
+ smu2 = numpy.linspace(-1, -0.01, ns2, dtype=float)
+ smu = numpy.array([(mu1, mu2) for mu1 in smu1 for mu2 in smu2])
+ return smu
+
+ def get_random_sample_of_mu(self, ns1: int = 1, ns2: int = 1):
+ smu = []
+ for i in range(ns1 * ns2):
+ smu1 = numpy.random.uniform(-1, -0.01)
+ smu2 = numpy.random.uniform(-1, -0.01)
+ smu.append((smu1, smu2))
+ smu = numpy.array(smu)
+ return smu
+
+ def get_bounds_on_space(self):
+ return [[min(self.x), max(self.x)], [min(self.y), max(self.y)]]
+
+ def get_bounds_on_parameter(self):
+ return [[-1, -0.01]] * 2
+
+ OneRealisation = FunctionH
+
+# ==============================================================================
+class LocalTest(unittest.TestCase):
+ @classmethod
+ def setUpClass(cls):
+ print('\nAUTODIAGNOSTIC\n==============\n')
+ print(" " + TwoDimensionalInverseDistanceCS2010().__doc__.strip())
+
+ def test001(self):
+ numpy.random.seed(123456789)
+ Equation = TwoDimensionalInverseDistanceCS2010()
+ for mu in Equation.get_sample_of_mu(5, 5):
+ solution = Equation.OneRealisation( mu )
+ # Nappe maximale au coin (0,0)
+ self.assertTrue(numpy.max(solution.flat) <= solution[0, 0])
+ # Nappe minimale au coin [-1,-1]
+ self.assertTrue(numpy.min(solution.flat) >= solution[-1, -1])
+
+ def tearDown(cls):
+ print("\n Tests OK\n")
+
+# ==============================================================================
+if __name__ == "__main__":
+ sys.stderr = sys.stdout
+ unittest.main(verbosity=0)
--- /dev/null
+# -*- coding: utf-8 -*-
+#
+# Copyright (C) 2008-2024 EDF R&D
+#
+# This library is free software; you can redistribute it and/or
+# modify it under the terms of the GNU Lesser General Public
+# License as published by the Free Software Foundation; either
+# version 2.1 of the License.
+#
+# This library is distributed in the hope that it will be useful,
+# but WITHOUT ANY WARRANTY; without even the implied warranty of
+# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
+# Lesser General Public License for more details.
+#
+# You should have received a copy of the GNU Lesser General Public
+# License along with this library; if not, write to the Free Software
+# Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
+#
+# See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+# Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
