From 1cadd6c5f233828ed3f98c566d9420c5c5018d6a Mon Sep 17 00:00:00 2001
From: Jean-Matthieu Etancelin <jean-matthieu.etancelin@imag.fr>
Date: Fri, 6 Jul 2012 14:16:32 +0000
Subject: [PATCH] Cleaning parmepy (in progress ...)
---
Examples/mainJM.py | 101 +++++----
HySoP/hysop/ParticularSolvers/basic.py | 173 ++++++++-------
HySoP/hysop/ParticularSolvers/solver.py | 9 +-
HySoP/hysop/domain/box.py | 3 +
HySoP/hysop/fields/__init__.py | 6 +
HySoP/hysop/fields/continuous.py | 23 +-
HySoP/hysop/fields/discrete.py | 52 ++++-
.../operator/{Advection.py => Transport.py} | 16 +-
HySoP/hysop/operator/advection_d.py | 198 +++++++++---------
HySoP/hysop/operator/continuous.py | 15 +-
HySoP/hysop/operator/discrete.py | 8 -
HySoP/hysop/problem/ContinuousProblem.py | 130 ------------
.../problem/ContinuousTransportProblem.py | 90 --------
HySoP/hysop/problem/__init__.py | 10 +-
HySoP/hysop/problem/continuous.py | 79 +++++++
.../{DiscreteProblem.py => discrete.py} | 0
HySoP/hysop/problem/problem.py | 90 ++++++++
17 files changed, 516 insertions(+), 487 deletions(-)
rename HySoP/hysop/operator/{Advection.py => Transport.py} (65%)
delete mode 100644 HySoP/hysop/problem/ContinuousProblem.py
delete mode 100644 HySoP/hysop/problem/ContinuousTransportProblem.py
create mode 100644 HySoP/hysop/problem/continuous.py
rename HySoP/hysop/problem/{DiscreteProblem.py => discrete.py} (100%)
create mode 100644 HySoP/hysop/problem/problem.py
diff --git a/Examples/mainJM.py b/Examples/mainJM.py
index 3e262b37b..31370b43b 100644
--- a/Examples/mainJM.py
+++ b/Examples/mainJM.py
@@ -1,6 +1,6 @@
# -*- coding: utf-8 -*-
import time
-import new_ParMePy as pp
+import parmepy as pp
def vitesse(x):
@@ -16,56 +16,73 @@ def scalaire(x):
def run():
# Parameter
- dim = 3
- nb = 128
- boxlength = [1., 1., 1.]
- boxMin = [0., 0., 0.]
- nbElem = [nb, nb, nb]
-
- timeStep = 0.02
- FinalTime = 1.
- outputFilePrefix = './res/RK2_'
- outputModulo = 0
+ # timeStep = 0.02
+ # FinalTime = 1.
+ # outputFilePrefix = './res/RK2_'
+ # outputModulo = 0
t0 = time.time()
- box = pp.Box(dimension=dim,
- length=boxlength,
- minPosition=boxMin)
- velo = pp.AnalyticalVariable(domain=box,
- dimension=dim,
- formula=vitesse, name='Velocity')
- scal = pp.AnalyticalVariable(domain=box,
- dimension=1,
- formula=scalaire, scalar=True, name='Scalar')
- advec = pp.Advection(velo, scal)
- cTspPb = pp.ContinuousTransportProblem(advection=advec)
-
- cTspPb.discretize(dSpec=[nbElem])
-
- InterpolationMethod = pp.Linear(grid=box.discreteDomain, gvalues=velo.discreteVariable)
- RemeshingMethod = pp.M6Prime(grid=box.discreteDomain, gvalues=scal.discreteVariable)
- ODESolver = pp.RK2(InterpolationMethod, box.discreteDomain.applyConditions)
- ## cpu
- #solver = pp.ParticularSolver(cTspPb, ODESolver, InterpolationMethod, RemeshingMethod)
- ## gpu
- solver = pp.GPUParticularSolver(cTspPb,
- ODESolver, InterpolationMethod, RemeshingMethod,
- device_type='gpu', src='./examples/mainJM_kernels.cl'
- )
- #solver = pp.GPUParticularSolver_GLRender(cTspPb, ODESolver, InterpolationMethod, RemeshingMethod, device_type='gpu', dt=timeStep)
- cTspPb.setSolver(solver)
- cTspPb.setPrinter(pp.Printer(frequency=outputModulo, outputPrefix=outputFilePrefix, problem=cTspPb))
+ ## Domain
+ box = pp.Box(3, length=[1., 1., 1.], origin=[0., 0., 0.])
+
+ ## Fields
+ scal = pp.ContinuousField(domain=box, name='Scalar')
+ velo = pp.ContinuousField(domain=box, name='Velocity', vector=True)
+
+ ## Operator
+ advec = pp.Transport(velo, scal)
+
+ ## Solver creation (discretisation of objects is done in solver initialisation)
+ topo3D = pp.CartesianTopology(domain=box, resolution=[4, 4, 4], dim=3)
+
+ ##Problem
+ pb = pp.Problem(topo3D, [advec])
+
+ ## Setting solver to Problem
+ pb.setSolver(pp.ParticularSolver)
+
+ # ## Setting io to Problem
+ # pb.setSolver(pp.Printer(...))
t1 = time.time()
- cTspPb.solve(T=FinalTime, dt=timeStep)
+ # ## Solve problem
+ # pb.solve(T=FinalTime, dt=timeStep)
- tf = time.time()
+ # tf = time.time()
- print "Total time : ", tf - t0, "sec (CPU)"
+ print pb
+ print "\n TODO : -terminer l'initialisation du solver (opérateurs) \n\t -Ajouter les opérateurs annexes (timer et io) au probleme \n\t -Ecrire méthode de résolution"
+
+ # print "Total time : ", tf - t0, "sec (CPU)"
print "Init time : ", t1 - t0, "sec (CPU)"
- print "Solving time : ", tf - t1, "sec (CPU)"
+ # print "Solving time : ", tf - t1, "sec (CPU)"
+
+
+ # InterpolationMethod = pp.Linear(grid=box.discreteDomain, gvalues=velo.discreteVariable)
+ # RemeshingMethod = pp.M6Prime(grid=box.discreteDomain, gvalues=scal.discreteVariable)
+ # ODESolver = pp.RK2(InterpolationMethod, box.discreteDomain.applyConditions)
+ # ## cpu
+ # #solver = pp.ParticularSolver(cTspPb, ODESolver, InterpolationMethod, RemeshingMethod)
+ # ## gpu
+ # solver = pp.GPUParticularSolver(cTspPb,
+ # ODESolver, InterpolationMethod, RemeshingMethod,
+ # device_type='gpu', src='./examples/mainJM_kernels.cl'
+ # )
+ # #solver = pp.GPUParticularSolver_GLRender(cTspPb, ODESolver, InterpolationMethod, RemeshingMethod, device_type='gpu', dt=timeStep)
+ # cTspPb.setSolver(solver)
+ # cTspPb.setPrinter(pp.Printer(frequency=outputModulo, outputPrefix=outputFilePrefix, problem=cTspPb))
+
+ # t1 = time.time()
+
+ # cTspPb.solve(T=FinalTime, dt=timeStep)
+
+ # tf = time.time()
+
+ # print "Total time : ", tf - t0, "sec (CPU)"
+ # print "Init time : ", t1 - t0, "sec (CPU)"
+ # print "Solving time : ", tf - t1, "sec (CPU)"
if __name__ == "__main__":
diff --git a/HySoP/hysop/ParticularSolvers/basic.py b/HySoP/hysop/ParticularSolvers/basic.py
index 714c18f9e..d6f7189b1 100644
--- a/HySoP/hysop/ParticularSolvers/basic.py
+++ b/HySoP/hysop/ParticularSolvers/basic.py
@@ -3,16 +3,8 @@
@package Utils
Particular solver description.
"""
-from ..Param import *
-from Solver import Solver
-from ..Domain.ParticleField import ParticleField
-from ..Variable.DiscreteVariable import DiscreteVariable
-from ..Operator.RemeshingDOp import RemeshingDOp
-from ..Operator.TagDOp import TagDOp
-from ..Operator.InterpolationDOp import InterpolationDOp
-from ..Operator.RemeshingDOp import RemeshingDOp
-
-from Lambda2 import Lambda2
+from solver import Solver
+from ..fields.continuous import ContinuousField
class ParticularSolver(Solver):
@@ -21,71 +13,75 @@ class ParticularSolver(Solver):
Link with differents numericals methods used.
"""
def __init__(self, problem,
- ODESolver,
- InterpolationMethod,
- RemeshingMethod,
+ ODESolver = None,
+ InterpolationMethod = None,
+ RemeshingMethod = None,
splittingMethod='strang'):
"""
Constructor.
Create a solver description.
"""
Solver.__init__(self, problem)
- self.discreteProblem = problem.discreteProblem
- ## ODE Solver method.
- self.ODESolver = ODESolver
- ## Interpolation method.
- self.InterpolationMethod = InterpolationMethod
- ## Remeshing method.
- self.RemeshingMethod = RemeshingMethod
- ## Splitting Method
- self.splittingMethod = splittingMethod
- self.dim = self.discreteProblem.domains[0].dimension
- self.grid = self.discreteProblem.domains[0]
- self.gvelo = self.discreteProblem.advecOp.velocity
- self.gscal = self.discreteProblem.advecOp.scalar
- self.parts = ParticleField(self.grid)
- self.ppos = DiscreteVariable(domain=self.parts, dimension=self.gvelo.dimension, name="Particle position")
- self.parts.setPositions(self.ppos)
- self.pvelo = DiscreteVariable(domain=self.parts, dimension=self.gvelo.dimension, name="Particle velocity")
- self.pscal = DiscreteVariable(domain=self.parts, dimension=self.gscal.dimension, scalar=True, name="Particle scalar")
- #self.pbloc = DiscreteVariable(domain=self.parts, dimension=1, integer=True, initialValue=lambda x: 0, scalar=True)
- #self.ptag = DiscreteVariable(domain=self.parts, dimension=1, integer=True, initialValue=lambda x: 0, scalar=True)
- # Interpolation : particles positions, grid velocity -> particles velocity
- interpol = InterpolationDOp(position=self.ppos,
- velocityField=self.discreteProblem.advecOp.velocity,
- resvelo=self.pvelo,
- method=self.InterpolationMethod)
- if isinstance(self.RemeshingMethod, Lambda2):
- # Tag particles : particles positions, particles velocities -> particles bloc, particles Tag
- tag = TagDOp(partPositions=self.ppos,
- partVelocities=self.pvelo,
- grid=self.discreteProblem.domains[0],
- resbloc=self.pbloc,
- restag=self.ptag)
- # Remeshing : particles positions, particles scalar, particles tag, particles bloc -> grid scalar
- self.remesh = RemeshingDOp(partPositions=self.ppos,
- partScalar=self.pscal,
- resscal=self.discreteProblem.advecOp.scalar,
- method=self.RemeshingMethod)
- # Advection : grid position, grid velocity, grid scalar -> particles position, particles scalar
- self.discreteProblem.advecOp.setResultVariable(respos=self.ppos, resscal=self.pscal)
- self.discreteProblem.advecOp.setMethod(self.ODESolver)
- if isinstance(self.RemeshingMethod, Lambda2):
- self.discreteProblem.addOperator([interpol, tag, self.remesh])
- else:
- self.discreteProblem.addOperator([self.remesh])
+ self.problem = problem
+ self.p_position = ContinuousField(domain=problem.domains[0], name='Particle Position')
+ self.p_scalar = ContinuousField(domain=problem.domains[0], name='Particle Scalar')
+ self.problem.addVariable([self.p_position, self.p_scalar])
+ # self.discreteProblem = problem.discreteProblem
+ # ## ODE Solver method.
+ # self.ODESolver = ODESolver
+ # ## Interpolation method.
+ # self.InterpolationMethod = InterpolationMethod
+ # ## Remeshing method.
+ # self.RemeshingMethod = RemeshingMethod
+ # ## Splitting Method
+ # self.splittingMethod = splittingMethod
+ # self.dim = self.discreteProblem.domains[0].dimension
+ # self.grid = self.discreteProblem.domains[0]
+ # self.gvelo = self.discreteProblem.advecOp.velocity
+ # self.gscal = self.discreteProblem.advecOp.scalar
+ # self.parts = ParticleField(self.grid)
+ # self.ppos = DiscreteVariable(domain=self.parts, dimension=self.gvelo.dimension, name="Particle position")
+ # self.parts.setPositions(self.ppos)
+ # self.pvelo = DiscreteVariable(domain=self.parts, dimension=self.gvelo.dimension, name="Particle velocity")
+ # self.pscal = DiscreteVariable(domain=self.parts, dimension=self.gscal.dimension, scalar=True, name="Particle scalar")
+ # #self.pbloc = DiscreteVariable(domain=self.parts, dimension=1, integer=True, initialValue=lambda x: 0, scalar=True)
+ # #self.ptag = DiscreteVariable(domain=self.parts, dimension=1, integer=True, initialValue=lambda x: 0, scalar=True)
+ # # Interpolation : particles positions, grid velocity -> particles velocity
+ # interpol = InterpolationDOp(position=self.ppos,
+ # velocityField=self.discreteProblem.advecOp.velocity,
+ # resvelo=self.pvelo,
+ # method=self.InterpolationMethod)
+ # if isinstance(self.RemeshingMethod, Lambda2):
+ # # Tag particles : particles positions, particles velocities -> particles bloc, particles Tag
+ # tag = TagDOp(partPositions=self.ppos,
+ # partVelocities=self.pvelo,
+ # grid=self.discreteProblem.domains[0],
+ # resbloc=self.pbloc,
+ # restag=self.ptag)
+ # # Remeshing : particles positions, particles scalar, particles tag, particles bloc -> grid scalar
+ # self.remesh = RemeshingDOp(partPositions=self.ppos,
+ # partScalar=self.pscal,
+ # resscal=self.discreteProblem.advecOp.scalar,
+ # method=self.RemeshingMethod)
+ # # Advection : grid position, grid velocity, grid scalar -> particles position, particles scalar
+ # self.discreteProblem.advecOp.setResultVariable(respos=self.ppos, resscal=self.pscal)
+ # self.discreteProblem.advecOp.setMethod(self.ODESolver)
+ # if isinstance(self.RemeshingMethod, Lambda2):
+ # self.discreteProblem.addOperator([interpol, tag, self.remesh])
+ # else:
+ # self.discreteProblem.addOperator([self.remesh])
- # Splitting Step TODO:Transformer le splitting en opérateur d'opérateur
- splitting = []
- if self.splittingMethod == 'strang':
- [splitting.append((i, 0.5)) for i in xrange(self.dim - 1)]
- splitting.append((self.dim - 1, 1.))
- [splitting.append((self.dim - 2 - i, 0.5)) for i in xrange(self.dim - 1)]
- else:
- raise ValueError("Unknown splitting method : ", self.splittingMethod)
- self.discreteProblem.operator_list = []
- [self.discreteProblem.operator_list.append((op, 1.0, 0)) for op in self.discreteProblem.operators if not op.is_splittable]
- [[self.discreteProblem.operator_list.append((op, i[1], i[0])) for op in self.discreteProblem.operators if op.is_splittable] for i in splitting]
+ # # Splitting Step TODO:Transformer le splitting en opérateur d'opérateur
+ # splitting = []
+ # if self.splittingMethod == 'strang':
+ # [splitting.append((i, 0.5)) for i in xrange(self.dim - 1)]
+ # splitting.append((self.dim - 1, 1.))
+ # [splitting.append((self.dim - 2 - i, 0.5)) for i in xrange(self.dim - 1)]
+ # else:
+ # raise ValueError("Unknown splitting method : ", self.splittingMethod)
+ # self.discreteProblem.operator_list = []
+ # [self.discreteProblem.operator_list.append((op, 1.0, 0)) for op in self.discreteProblem.operators if not op.is_splittable]
+ # [[self.discreteProblem.operator_list.append((op, i[1], i[0])) for op in self.discreteProblem.operators if op.is_splittable] for i in splitting]
def initialize(self):
"""
@@ -93,32 +89,29 @@ class ParticularSolver(Solver):
@param discreteProblem : Problem to initialize.
"""
- ### Create grid.points TODO: A enlever
- self.grid.points = np.empty_like(self.gvelo.values, dtype=dtype_real)
- if self.dim == 1:
- for x in xrange(self.grid.elementNumber[0]):
- self.grid.points[x] = self.grid.axes[0][x]
- if self.dim == 2:
- for x in xrange(self.grid.elementNumber[0]):
- for y in xrange(self.grid.elementNumber[1]):
- self.grid.points[x, y][0] = self.grid.axes[0][x]
- self.grid.points[x, y][1] = self.grid.axes[1][y]
- if self.dim == 3:
- for x in xrange(self.grid.elementNumber[0]):
- for y in xrange(self.grid.elementNumber[1]):
- for z in xrange(self.grid.elementNumber[2]):
- self.grid.points[x, y, z][0] = self.grid.axes[0][x]
- self.grid.points[x, y, z][1] = self.grid.axes[1][y]
- self.grid.points[x, y, z][2] = self.grid.axes[2][z]
- self.gvelo.init()
- self.gscal.init()
+ # ### Create grid.points TODO: A enlever
+ # self.grid.points = np.empty_like(self.gvelo.values, dtype=dtype_real)
+ # if self.dim == 1:
+ # for x in xrange(self.grid.elementNumber[0]):
+ # self.grid.points[x] = self.grid.axes[0][x]
+ # if self.dim == 2:
+ # for x in xrange(self.grid.elementNumber[0]):
+ # for y in xrange(self.grid.elementNumber[1]):
+ # self.grid.points[x, y][0] = self.grid.axes[0][x]
+ # self.grid.points[x, y][1] = self.grid.axes[1][y]
+ # if self.dim == 3:
+ # for x in xrange(self.grid.elementNumber[0]):
+ # for y in xrange(self.grid.elementNumber[1]):
+ # for z in xrange(self.grid.elementNumber[2]):
+ # self.grid.points[x, y, z][0] = self.grid.axes[0][x]
+ # self.grid.points[x, y, z][1] = self.grid.axes[1][y]
+ # self.grid.points[x, y, z][2] = self.grid.axes[2][z]
+ # self.gvelo.init()
+ # self.gscal.init()
def __str__(self):
"""ToString method"""
s = " Particular solver "
- s += "\n - ODESolver : " + str(self.ODESolver)
- s += "\n - Interpolation : " + str(self.InterpolationMethod)
- s += "\n - Remeshing : " + str(self.RemeshingMethod)
return s
if __name__ == "__main__":
diff --git a/HySoP/hysop/ParticularSolvers/solver.py b/HySoP/hysop/ParticularSolvers/solver.py
index d0a2ddbd1..0f133c266 100644
--- a/HySoP/hysop/ParticularSolvers/solver.py
+++ b/HySoP/hysop/ParticularSolvers/solver.py
@@ -1,21 +1,26 @@
# -*- coding: utf-8 -*-
"""
-@package Utils
+@package solver
Solver interface.
"""
-from ..Param import *
+from abc import ABCMeta, abstractmethod
class Solver:
"""
Solver interface.
"""
+
+ __metaclass__ = ABCMeta
+
+ @abstractmethod
def __init__(self, problem):
"""
Constructor of an interface. Do nothing.
"""
pass
+ @abstractmethod
def initialize(self):
"""
Solver initialisation for a given DiscreteProblem.DiscreteProblem.
diff --git a/HySoP/hysop/domain/box.py b/HySoP/hysop/domain/box.py
index 6ba485a0d..ed0528797 100644
--- a/HySoP/hysop/domain/box.py
+++ b/HySoP/hysop/domain/box.py
@@ -33,6 +33,7 @@ class Box(Domain):
# Boundary conditions type :
self.boundaries = np.zeros((self.dimension))
self.boundaries[:] = PERIODIC
+ self.discreteDomain = None
def discretize(self, resolution):
"""Box discretization method.
@@ -47,6 +48,8 @@ class Box(Domain):
""" Doc display. """
s = str(self.dimension) + "D box (parallelepipedic or rectangular) domain : \n"
s += " origin : " + str(self.origin) + ", maxPosition :" + str(self.max) + ", lengthes :" + str(self.length) + "."
+ if self.discreteDomain is not None:
+ s += "Discretised on : \n" + str(self.discreteDomain)
return s
if __name__ == "__main__":
diff --git a/HySoP/hysop/fields/__init__.py b/HySoP/hysop/fields/__init__.py
index e69de29bb..cee2ba565 100644
--- a/HySoP/hysop/fields/__init__.py
+++ b/HySoP/hysop/fields/__init__.py
@@ -0,0 +1,6 @@
+"""
+@package parmepy.fields
+
+Everything concerning Fields.
+
+"""
diff --git a/HySoP/hysop/fields/continuous.py b/HySoP/hysop/fields/continuous.py
index 09d8b5909..284d4d404 100644
--- a/HySoP/hysop/fields/continuous.py
+++ b/HySoP/hysop/fields/continuous.py
@@ -4,12 +4,13 @@
Continuous variable description
"""
from .discrete import ScalarField
+from .discrete import VectorField
class ContinuousField(object):
""" A variable defined on a physical domain """
- def __init__(self, domain, name="?"):
+ def __init__(self, domain, name="?", vector=False):
"""
Constructor.
Create a ContinuousField.
@@ -29,6 +30,8 @@ class ContinuousField(object):
self.__fieldId = -1
## list of the various discretizations of this field
self.discreteField = []
+ ## is field is a vector field
+ self.vector = vector
def discretize(self, topology=None):
"""
@@ -36,18 +39,28 @@ class ContinuousField(object):
"""
self.__fieldId += 1
- self.discreteField.append(ScalarField(self, topology, name=self.name + "_D_" + str(self.__fieldId),
- idFromParent=self.__fieldId))
+ print "Discretisation de ", self.name,
+ if self.vector:
+ print "Vector Field"
+ self.discreteField.append(VectorField(self, topology,
+ name=self.name + "_D_" + str(self.__fieldId),
+ idFromParent=self.__fieldId))
+ else:
+ print "Scalar Field"
+ self.discreteField.append(ScalarField(self, topology,
+ name=self.name + "_D_" + str(self.__fieldId),
+ idFromParent=self.__fieldId))
return self.discreteField[self.__fieldId], self.__fieldId
def __str__(self):
""" Field info display """
- s = self.name + ' , ' + str(self.dimension) + 'D field with the following discretizations:\n'
+ s = self.name + ' , ' + str(self.dimension) + 'D field '
if len(self.discreteField) != 0:
+ s += "with the following discretizations:\n"
for i in range(len(self.discreteField)):
s += self.discreteField[i].__str__()
else:
- s += "None."
+ s += ". Not discretised\n"
return s
diff --git a/HySoP/hysop/fields/discrete.py b/HySoP/hysop/fields/discrete.py
index 408351a5e..8f05381c2 100644
--- a/HySoP/hysop/fields/discrete.py
+++ b/HySoP/hysop/fields/discrete.py
@@ -83,33 +83,63 @@ class ScalarField(object):
class VectorField(object):
"""
A vector field, defined on each grid of each subdomain of the given topology.
-
+ Vector field is composed of several ScalarFields
"""
- def __init__(self, topology, name=""):
- self.topology = topology
+ def __init__(self, parent, topology=None, name="?", idFromParent=None):
+ if(topology is not None):
+ self.topology = topology
+ else:
+ raise NotImplementedError()
+
self.name = name
self.dimension = topology.domain.dimension
resolution = self.topology.mesh.resolution
- self.data = np.zeros((resolution), dtype=PARMES_REAL, order=ORDER)
+ ## \todo is numpy array zeros or empty???
+ self.data = [np.zeros((resolution), dtype=PARMES_REAL, order=ORDER) for d in xrange(self.dimension)]
+ self.__parentField = parent
+ self.__idFromParent = idFromParent
def __str__(self):
""" Display class information """
- s = self.name + ": " + str(self.dimension) + "D scalar field of shape " + str(self.data.shape)
+ s = '[' + str(self.topology.rank) + '] '
+ s += " id " + str(self.__idFromParent) + ", "
+ s += self.name + " " + str(self.dimension) + 'D discrete vector field with resolution '
+ s += str([data.shape for data in self.data]) + "."
return s + "\n"
def __getitem__(self, i):
- """ Access to the content of the field """
- return self.data.__getitem__(i)
+ """ Access to the content of the field
+ Usage (3D): A = VectorField(...), We have 3 components len(a.data) == 3.
+ Following instructions access to index 2,1,1 of y component:
+ A[2,1,1,1]
+ A[1][2,1,1]
+ Access to whole vector of index 2,1,1: A[2,1,1]
+ """
+ try:
+ if len(i) == len(self.data):
+ return [data.__getitem__(i) for data in self.data]
+ else:
+ return self.data[i[-1]].__getitem__(tuple(i[0:-1]))
+ except (TypeError):
+ return self.data[i]
def __setitem__(self, i, value):
"""
- Set the content of the field component at position i
+ Set the content of the vector field component at position i
Usage :
- A = ScalarField(topo)
- A[2,1,1] = 12.0 # Calls A.data[2,1,1] = 12.0
+ A = VectorField(topo)
+ A[2,1,1] = 12.0 # Calls A.data[d][2,1,1] = 12.0 for all components d.
+ A[2,1,1] = [12.0, 13.0, 14.0] # Calls A.data[0][2,1,1] = 12.0, A.data[1][2,1,1] = 13.0 and A.data[2][2,1,1] = 14.0
+ A[1][2,1,1] = 13.0 # Calls A.data[1][2,1,1] = 12.0
"""
- self.data.__setitem__(i, value)
+ if len(i) == len(self.data):
+ try:
+ [data.__setitem__(i, v) for data, v in zip(self.data, value)]
+ except (TypeError):
+ [data.__setitem__(i, value) for data in self.data]
+ else:
+ self.data[i[-1]].__getitem__(tuple(i[0:-1]))
def initialize(self, formula):
"""
diff --git a/HySoP/hysop/operator/Advection.py b/HySoP/hysop/operator/Transport.py
similarity index 65%
rename from HySoP/hysop/operator/Advection.py
rename to HySoP/hysop/operator/Transport.py
index 144b51902..00541de94 100644
--- a/HySoP/hysop/operator/Advection.py
+++ b/HySoP/hysop/operator/Transport.py
@@ -5,10 +5,10 @@
Advection operator representation
"""
from .continuous import ContinuousOperator
-from .advection_d import AdvectionDOp
+from .advection_d import Advection_P
-class Advection(ContinuousOperator):
+class Transport(ContinuousOperator):
"""
Advection operator representation
"""
@@ -26,8 +26,9 @@ class Advection(ContinuousOperator):
self.velocity = velocity
## advected scalar variable
self.scalar = scalar
+ self.addVariable([velocity, scalar])
- def discretize(self, spec=None):
+ def discretize(self, idVelocityD=0, idScalarD=0, result_position=None, result_scalar=None):
"""
Advection operator discretization method.
Create an AdvectionDOp.AdvectionDOp from given specifications.
@@ -35,11 +36,16 @@ class Advection(ContinuousOperator):
@param spec : discretization specifications, not used.
"""
if self.discreteOperator is None:
- self.discreteOperator = AdvectionDOp(self)
+ self.discreteOperator = Advection_P(self, idVelocityD, idScalarD, result_position, result_scalar)
def __str__(self):
"""ToString method"""
- s = " Advection operator (ContinuousOperator)"
+ s = "Advection operator (ContinuousOperator)"
+ if self.discreteOperator is not None:
+ s += "with the following discretization:\n"
+ s += str(self.discreteOperator)
+ else:
+ s += ". Not discretised"
return s + "\n"
if __name__ == "__main__":
diff --git a/HySoP/hysop/operator/advection_d.py b/HySoP/hysop/operator/advection_d.py
index 8aaab8ea2..94084acd9 100644
--- a/HySoP/hysop/operator/advection_d.py
+++ b/HySoP/hysop/operator/advection_d.py
@@ -8,122 +8,124 @@ import pyopencl as cl
import time
-class AdvectionDOp(DiscreteOperator):
+class Advection_P(DiscreteOperator):
"""
- Advection operator representation.
+ Particle avection operator representation.
DiscreteOperator.DiscreteOperator specialization.
"""
- def __init__(self, advec):
+ def __init__(self, advec, idVelocityD=0, idScalarD=0, result_position=None, result_scalar=None):
"""
Constructor.
- Create a Transport operator on a given continuous domain.
+ Create a Advection operator on a given continuous domain.
Work on a given scalar at a given velocity to produce scalar distribution at new positions.
@param advec : Advection operator.
"""
DiscreteOperator.__init__(self)
- ## Velocity.
- self.velocity = advec.velocity.discreteVariable
+ ## Velocity. \todo utiliser un vectorField
+ self.velocity = advec.velocity.discreteField[idVelocityD]
## Transported scalar.
- self.scalar = advec.scalar.discreteVariable
- self.addVariable([self.velocity, self.scalar])
- self.gpu_kernel_name = "advection"
- self.compute_time = [0., 0., 0.]
- self.submit_time = 0.
- self.queued_time = 0.
- self.total_flop = 0
- self.total_bytes_accessed = 0
-
- def setResultVariable(self, respos, resscal):
- """
- Set results variables
-
- @param respos DiscreteVariable.DiscreteVariable : new positions
- @param resscal DiscreteVariable.DiscreteVariable : scalar values
- """
- ## Positions after advection
- self.resultPosition = respos
- ## Scalar value on new positions
- self.resultScalar = resscal
+ self.scalar = advec.scalar.discreteField[idScalarD]
+ ## Result position
+ self.res_position = result_position
+ ## Result scalar
+ self.res_scalar = result_scalar
+ ## input fields
+ self.input = [self.velocity, self.scalar]
+ ## output fields
+ self.output = [self.res_position, self.res_scalar]
+ # self.gpu_kernel_name = "advection"
+ # self.compute_time = [0., 0., 0.]
+ # self.submit_time = 0.
+ # self.queued_time = 0.
+ # self.total_flop = 0
+ # self.total_bytes_accessed = 0
def applyOperator(self, t, dt, splittingDirection):
- """
- Apply advection operator.
+ raise NotImplementedMethod("Not yet implemented")
+ # """
+ # Apply advection operator.
- @param t : current time.
- @param dt : time step for advection.
- @param splittingDirection : direction to advect.
- """
- if self.gpu_kernel is None:
- c_time = time.time()
- ## ------------------- NumPy Optimisation
- # Particle update
- self.resultScalar.values[...] = self.scalar.values[...]
- self.resultPosition.values[...] = self.velocity.domain.points[...]
- # Advection
- self.numMethod.integrate(self.velocity.domain.points, self.velocity.values, self.resultPosition.values, t, dt, splittingDirection)
- #print self.resultPosition.values
- ## -------------------
- self.compute_time[splittingDirection] += (time.time() - c_time)
- self.total_flop += self.numMethod.nb_flop * np.prod(self.resultPosition.values.shape[0:-1])
- else:
- l = list(self.gpu_shape)
- nb_ligne = l.pop(splittingDirection)
- local_wg = list(l)
- if len(l) == 2:
- local_wg[0] = 1
- local_wg = tuple(local_wg)
- else:
- local_wg = None
- if debug == 0:
- evt_copy = cl.enqueue_copy(self.gpu_queue, self.resultScalar.gpu_mem_object, self.scalar.gpu_mem_object)
- evt = self.gpu_kernel.advection(self.gpu_queue, tuple(l), local_wg,
- self.velocity.gpu_mem_object,
- self.scalar.gpu_mem_object,
- self.resultPosition.gpu_mem_object,
- self.resultScalar.gpu_mem_object,
- dtype_real_gpu(t), dtype_real_gpu(dt), dtype_integer(splittingDirection), dtype_integer(nb_ligne),
- self.velocity.domain.min.astype(dtype_real_gpu),
- self.velocity.domain.length.astype(dtype_real_gpu),
- self.velocity.domain.elementNumber.astype(dtype_integer),
- self.velocity.domain.elementSize.astype(dtype_real_gpu)
- )
- else:
- evt = self.gpu_kernel.advection(self.gpu_queue, tuple(l), None,
- self.velocity.gpu_mem_object,
- self.scalar.gpu_mem_object,
- self.resultPosition.gpu_mem_object,
- self.resultScalar.gpu_mem_object,
- self.debug_float_buffer, self.debug_integer_buffer,
- dtype_real_gpu(t), dtype_real_gpu(dt), dtype_integer(splittingDirection), dtype_integer(nb_ligne),
- self.velocity.domain.min.astype(dtype_real_gpu),
- self.velocity.domain.length.astype(dtype_real_gpu),
- self.velocity.domain.elementNumber.astype(dtype_integer),
- self.velocity.domain.elementSize.astype(dtype_real_gpu)
- )
- # self.gpu_queue.finish()
- # cl.enqueue_copy(self.gpu_queue, self.resultPosition.values, self.resultPosition.gpu_mem_object)
- # print self.resultPosition.values
- self.gpu_queue.finish()
- if debug > 0:
- cl.enqueue_copy(self.gpu_queue, self.debug_float, self.debug_float_buffer)
- cl.enqueue_copy(self.gpu_queue, self.debug_integer, self.debug_integer_buffer)
- print " :: DEBUG :: ", self.debug_float
- print " :: DEBUG :: ", self.debug_integer
- self.queued_time += (evt.profile.submit - evt.profile.queued) * 1e-9
- self.submit_time += (evt.profile.start - evt.profile.submit) * 1e-9
- self.compute_time[splittingDirection] += (evt.profile.end - evt.profile.start) * 1e-9
- self.queued_time += (evt_copy.profile.submit - evt_copy.profile.queued) * 1e-9
- self.submit_time += (evt_copy.profile.start - evt_copy.profile.submit) * 1e-9
- self.compute_time[splittingDirection] += (evt_copy.profile.end - evt_copy.profile.start) * 1e-9
- self.total_flop += 20 * nb_ligne * np.prod(l) # 20 flop par particule d'une ligne
- self.total_bytes_accessed += 2 * nb_ligne * np.prod(l) * 4 # 2 floats lus + 1+DIM float écrits
- #print "Advection:", self.compute_time * 1e-9
+ # @param t : current time.
+ # @param dt : time step for advection.
+ # @param splittingDirection : direction to advect.
+ # """
+ # if self.gpu_kernel is None:
+ # c_time = time.time()
+ # ## ------------------- NumPy Optimisation
+ # # Particle update
+ # self.resultScalar.values[...] = self.scalar.values[...]
+ # self.resultPosition.values[...] = self.velocity.domain.points[...]
+ # # Advection
+ # self.numMethod.integrate(self.velocity.domain.points, self.velocity.values, self.resultPosition.values, t, dt, splittingDirection)
+ # #print self.resultPosition.values
+ # ## -------------------
+ # self.compute_time[splittingDirection] += (time.time() - c_time)
+ # self.total_flop += self.numMethod.nb_flop * np.prod(self.resultPosition.values.shape[0:-1])
+ # else:
+ # l = list(self.gpu_shape)
+ # nb_ligne = l.pop(splittingDirection)
+ # local_wg = list(l)
+ # if len(l) == 2:
+ # local_wg[0] = 1
+ # local_wg = tuple(local_wg)
+ # else:
+ # local_wg = None
+ # if debug == 0:
+ # evt_copy = cl.enqueue_copy(self.gpu_queue, self.resultScalar.gpu_mem_object, self.scalar.gpu_mem_object)
+ # evt = self.gpu_kernel.advection(self.gpu_queue, tuple(l), local_wg,
+ # self.velocity.gpu_mem_object,
+ # self.scalar.gpu_mem_object,
+ # self.resultPosition.gpu_mem_object,
+ # self.resultScalar.gpu_mem_object,
+ # dtype_real_gpu(t), dtype_real_gpu(dt), dtype_integer(splittingDirection), dtype_integer(nb_ligne),
+ # self.velocity.domain.min.astype(dtype_real_gpu),
+ # self.velocity.domain.length.astype(dtype_real_gpu),
+ # self.velocity.domain.elementNumber.astype(dtype_integer),
+ # self.velocity.domain.elementSize.astype(dtype_real_gpu)
+ # )
+ # else:
+ # evt = self.gpu_kernel.advection(self.gpu_queue, tuple(l), None,
+ # self.velocity.gpu_mem_object,
+ # self.scalar.gpu_mem_object,
+ # self.resultPosition.gpu_mem_object,
+ # self.resultScalar.gpu_mem_object,
+ # self.debug_float_buffer, self.debug_integer_buffer,
+ # dtype_real_gpu(t), dtype_real_gpu(dt), dtype_integer(splittingDirection), dtype_integer(nb_ligne),
+ # self.velocity.domain.min.astype(dtype_real_gpu),
+ # self.velocity.domain.length.astype(dtype_real_gpu),
+ # self.velocity.domain.elementNumber.astype(dtype_integer),
+ # self.velocity.domain.elementSize.astype(dtype_real_gpu)
+ # )
+ # # self.gpu_queue.finish()
+ # # cl.enqueue_copy(self.gpu_queue, self.resultPosition.values, self.resultPosition.gpu_mem_object)
+ # # print self.resultPosition.values
+ # self.gpu_queue.finish()
+ # if debug > 0:
+ # cl.enqueue_copy(self.gpu_queue, self.debug_float, self.debug_float_buffer)
+ # cl.enqueue_copy(self.gpu_queue, self.debug_integer, self.debug_integer_buffer)
+ # print " :: DEBUG :: ", self.debug_float
+ # print " :: DEBUG :: ", self.debug_integer
+ # self.queued_time += (evt.profile.submit - evt.profile.queued) * 1e-9
+ # self.submit_time += (evt.profile.start - evt.profile.submit) * 1e-9
+ # self.compute_time[splittingDirection] += (evt.profile.end - evt.profile.start) * 1e-9
+ # self.queued_time += (evt_copy.profile.submit - evt_copy.profile.queued) * 1e-9
+ # self.submit_time += (evt_copy.profile.start - evt_copy.profile.submit) * 1e-9
+ # self.compute_time[splittingDirection] += (evt_copy.profile.end - evt_copy.profile.start) * 1e-9
+ # self.total_flop += 20 * nb_ligne * np.prod(l) # 20 flop par particule d'une ligne
+ # self.total_bytes_accessed += 2 * nb_ligne * np.prod(l) * 4 # 2 floats lus + 1+DIM float écrits
+ # #print "Advection:", self.compute_time * 1e-9
def __str__(self):
"""ToString method"""
- return "AdvectionDOp (DiscreteOperator)"
+ s = "Advection_P (DiscreteOperator). Applying on following fields : \n Input fields : \n"
+ for f in self.input:
+ s += str(f) + "\n"
+ s += "Output fields : \n"
+ for f in self.output:
+ s += str(f) + "\n"
+ return s
if __name__ == "__main__":
print __doc__
diff --git a/HySoP/hysop/operator/continuous.py b/HySoP/hysop/operator/continuous.py
index 6241f5617..ed2b750a0 100644
--- a/HySoP/hysop/operator/continuous.py
+++ b/HySoP/hysop/operator/continuous.py
@@ -20,11 +20,22 @@ class ContinuousOperator:
Constructor.
Create a ContinuousOperator.
"""
- ## Application domains of the variables.
- self.domains = []
+ ## Variables
+ self.variables = []
## Operator discretization.
self.discreteOperator = None
+ def addVariable(self, cVariable):
+ """
+ Add an continuous variables to the operator.
+ Also add variables' domains to the operator.
+
+ @param cVariable ContinuousVariable.ContinuousVariable : variables to add.
+ """
+ for v in cVariable:
+ if self.variables.count(v) == 0:
+ self.variables.append(v)
+
@abstractmethod
def discretize(self, spec=None):
"""
diff --git a/HySoP/hysop/operator/discrete.py b/HySoP/hysop/operator/discrete.py
index d816c828b..2578a9755 100644
--- a/HySoP/hysop/operator/discrete.py
+++ b/HySoP/hysop/operator/discrete.py
@@ -35,14 +35,6 @@ class DiscreteOperator:
"""
self.numMethod = method
- @abstractmethod
- def setResultVariable(self):
- """
- Abstract method, apply operaton on a variable.
- Must be implemented by sub-class.
- """
- raise NotImplementedError("Need to override method in a subclass of " + providedClass)
-
@abstractmethod
def applyOperator(self):
"""
diff --git a/HySoP/hysop/problem/ContinuousProblem.py b/HySoP/hysop/problem/ContinuousProblem.py
deleted file mode 100644
index cce6b951c..000000000
--- a/HySoP/hysop/problem/ContinuousProblem.py
+++ /dev/null
@@ -1,130 +0,0 @@
-# -*- coding: utf-8 -*-
-"""
-@package Problem
-Problem representation
-"""
-from ..Param import *
-
-
-class ContinuousProblem:
- """
- Abstract description of continuous problem.
- """
-
- def __init__(self):
- """
- Constructor.
- Create an empty problem.
- """
- ## Domains of the problem.
- self.domains = []
- ## Variables of the problem.
- self.variables = []
- ## Operators of the problem.
- self.operators = []
- ## Discretization of the problem.
- self.discreteProblem = None
-
- def addDomain(self, cDomain):
- """
- Add a list of ContinuousDomain.ContinuousDomain to the problem.
-
- @param cDomain ContinuousDomain.ContinuousDomain : list of domains to add.
- """
- for d in cDomain:
- if self.domains.count(d) == 0:
- self.domains.append(d)
-
- def addVariable(self, cVariable):
- """
- Add a list of ContinuousVariable.ContinuousVariable to the problem.
- Also add variables' domains to the problem.
-
- @param cVariable ContinuousVariable.ContinuousVariable : list of variables to add.
- """
- for v in cVariable:
- if self.variables.count(v) == 0:
- self.variables.append(v)
- if self.domains.count(v.domain) == 0:
- self.domains.append(v.domain)
-
- def addOperator(self, cOperator, index=None):
- """
- Add a list of ContinuousOperator.ContinuousOperator to the problem.
- Also add operators' variables and variables' domains to the problem.
-
- @param cOperator ContinuousOperator.ContinuousOperator : list of operators to add.
- """
- if isinstance(cOperator, list):
- for i,o in enumerate(cOperator):
- if self.operators.count(o) == 0:
- try:
- self.operators.insert(index[i],o)
- except TypeError:
- self.operators.append(o)
- for v in o.variables:
- if self.variables.count(v) == 0:
- self.variables.append(v)
- if self.domains.count(v.domain) == 0:
- self.domains.append(v.domain)
- else:
- if self.operators.count(cOperator) == 0:
- try:
- self.operators.insert(index, cOperator)
- except TypeError:
- self.operators.append(cOperator)
- for v in cOperator.variables:
- if self.variables.count(v) == 0:
- self.variables.append(v)
- if self.domains.count(v.domain) == 0:
- self.domains.append(v.domain)
-
- def discretize(self, dSpec=None, vSpec=None, oSpec=None):
- """
- Abstract method.
- Must be implemented by sub-class.
-
- @param dSpec list : ContinuousDomain.ContinuousDomain discretization specifications
- @param vSpec list : ContinuousVariable.ContinuousVariable discretization specifications
- @param oSpec list : ContinuousOperator.ContinuousOperator discretization specifications
- """
- raise NotImplementedError("Need to override method in a subclass of " + providedClass)
-
- def setSolver(self, solver):
- """
- Set solver for the discrete problem.
- It must have a discrete Problem.
-
- @param solver : Solver.Solver to use.
- """
- if self.discreteProblem == None:
- raise NotImplementedError("Cannot set solver on a non discretized problem " + providedClass)
- self.discreteProblem.setSolver(solver)
-
- def setPrinter(self, printer):
- """
- Set result printer.
- It must have a discrete Problem.
-
- @param printer : Printer.Printer to use.
- """
- if self.discreteProblem == None:
- raise NotImplementedError("Cannot set printer on a non discretized problem " + providedClass)
- self.discreteProblem.setPrinter(printer)
-
- def solve(self, T, dt):
- """
- Solving the problem.
- It must have a discrete Problem to solve.
-
- @param T : Simulation final time
- @param dt : Simulation time step
- """
- if self.discreteProblem == None:
- raise NotImplementedError("Cannot solve a non discretized problem " + providedClass)
- self.discreteProblem.solve(T, dt)
-
-if __name__ == "__main__":
- print __doc__
- print "- Provided class : ContinuousProblem"
- print ContinuousProblem.__doc__
diff --git a/HySoP/hysop/problem/ContinuousTransportProblem.py b/HySoP/hysop/problem/ContinuousTransportProblem.py
deleted file mode 100644
index 969a98eed..000000000
--- a/HySoP/hysop/problem/ContinuousTransportProblem.py
+++ /dev/null
@@ -1,90 +0,0 @@
-# -*- coding: utf-8 -*-
-"""
-@package Problem
-Problem representation
-"""
-from ..Param import *
-from ContinuousProblem import ContinuousProblem
-from DiscreteTransportProblem import DiscreteTransportProblem
-
-
-class ContinuousTransportProblem(ContinuousProblem):
- """
- Continuous Transport problem representation.
- """
-
- def __init__(self, advection, domains=None, variables=None, operators=None):
- """
- Constructor.
- Create a transport problem instance.
-
- @param advection : Advection operator for transport problem.
- @param domains : Additionals domains.
- @param variables : Additionnals variables.
- @param operators : Additionnals operators.
- """
- ContinuousProblem.__init__(self)
- ## Advection operator for the transport problem
- self.advecOp = advection
- self.addOperator([self.advecOp])
- if domains != None:
- self.addDomain(domains)
- if variables != None:
- self.addVariable(variables)
- if operators != None:
- self.addOperator(operators)
-
- def discretize(self, dSpec=None, vSpec=None, oSpec=None):
- """
- Transport operator discretization.
- Create a DiscreteTransportProblem from given discretizations specifications.
-
- @param dSpec list : ContinuousDomain.ContinuousDomain discretization specifications, element number for each direction.
- @param vSpec list : ContinuousVariable.ContinuousVariable discretization specifications, not used.
- @param oSpec list : ContinuousOperator.ContinuousOperator discretization specifications, not used.
- """
- if dSpec == None:
- dSpec = [None for d in self.domains]
- if vSpec == None:
- vSpec = [None for d in self.variables]
- if oSpec == None:
- oSpec = [None for d in self.operators]
- # Discretizations
- for i in xrange(len(self.domains)):
- self.domains[i].discretize(dSpec[i])
- for i in xrange(len(self.variables)):
- self.variables[i].discretize(vSpec[i])
- for i in xrange(len(self.operators)):
- self.operators[i].discretize(oSpec[i])
- self.discreteProblem = DiscreteTransportProblem(advection=self.advecOp.discreteOperator)
- self.discreteProblem.addDomain([d.discreteDomain for d in self.domains])
- self.discreteProblem.addVariable([v.discreteVariable for v in self.variables])
- self.discreteProblem.addOperator([o.discreteOperator for o in self.operators], [self.operators.index(o) for o in self.operators])
-
- def writeToFile(self, fileStr):
- """
- Write the data to a file.
-
- @param fileStr String : file name
- """
- self.discreteProblem.writeToFile(fileStr)
-
- def __str__(self):
- """ToString method"""
- s = "ContinuousTransportProblem (ContinuousProblem)\n"
- s += " {0} domains, {1} variables and {2} operators.\n".format(len(self.domains), len(self.variables), len(self.operators))
- s += " Domains : \n"
- for d in self.domains:
- s += str(d)
- s += "\n Variables : \n"
- for v in self.variables:
- s += str(v)
- s += "\n Operators : \n"
- for o in self.operators:
- s += str(o)
- return s
-
-if __name__ == "__main__":
- print __doc__
- print "- Provided class : ContinuousTransportProblem"
- print ContinuousTransportProblem.__doc__
diff --git a/HySoP/hysop/problem/__init__.py b/HySoP/hysop/problem/__init__.py
index e6d9b54b3..1553dd6dc 100644
--- a/HySoP/hysop/problem/__init__.py
+++ b/HySoP/hysop/problem/__init__.py
@@ -1,4 +1,6 @@
-from ContinuousProblem import ContinuousProblem
-from ContinuousTransportProblem import ContinuousTransportProblem
-from DiscreteProblem import DiscreteProblem
-from DiscreteTransportProblem import DiscreteTransportProblem
+"""
+@package parmepy.problem
+
+Everything concerning Problems.
+
+"""
diff --git a/HySoP/hysop/problem/continuous.py b/HySoP/hysop/problem/continuous.py
new file mode 100644
index 000000000..125e63005
--- /dev/null
+++ b/HySoP/hysop/problem/continuous.py
@@ -0,0 +1,79 @@
+# -*- coding: utf-8 -*-
+"""
+@package Problem
+Problem representation
+"""
+from abc import ABCMeta, abstractmethod
+
+
+class ContinuousProblem:
+ """
+ Abstract description of continuous problem.
+ """
+
+ __metaclass__ = ABCMeta
+
+ @abstractmethod
+ def __init__(self):
+ """
+ Constructor.
+ Create an empty problem.
+ """
+ ## Operators of the problem.
+ self.operators = []
+ ## Discretization of the problem.
+ self.discreteProblem = None
+
+ @abstractmethod
+ def discretize(self, spec=None):
+ """
+ Abstract method.
+ Must be implemented by sub-class.
+
+ @param dSpec list : ContinuousDomain.ContinuousDomain discretization specifications
+ @param vSpec list : ContinuousVariable.ContinuousVariable discretization specifications
+ @param oSpec list : ContinuousOperator.ContinuousOperator discretization specifications
+ """
+ raise NotImplementedError("Need to override method in a subclass")
+
+ @abstractmethod
+ def setSolver(self, solver):
+ """
+ Set solver for the discrete problem.
+ It must have a discrete Problem.
+
+ @param solver : Solver.Solver to use.
+ """
+ if self.discreteProblem == None:
+ raise NotImplementedError("Cannot set solver on a non discretized problem ")
+ self.discreteProblem.setSolver(solver)
+
+ @abstractmethod
+ def setPrinter(self, printer):
+ """
+ Set result printer.
+ It must have a discrete Problem.
+
+ @param printer : Printer.Printer to use.
+ """
+ if self.discreteProblem == None:
+ raise NotImplementedError("Cannot set printer on a non discretized problem ")
+ self.discreteProblem.setPrinter(printer)
+
+ @abstractmethod
+ def solve(self, T, dt):
+ """
+ Solving the problem.
+ It must have a discrete Problem to solve.
+
+ @param T : Simulation final time
+ @param dt : Simulation time step
+ """
+ if self.discreteProblem == None:
+ raise NotImplementedError("Cannot solve a non discretized problem ")
+ self.discreteProblem.solve(T, dt)
+
+if __name__ == "__main__":
+ print __doc__
+ print "- Provided class : ContinuousProblem"
+ print ContinuousProblem.__doc__
diff --git a/HySoP/hysop/problem/DiscreteProblem.py b/HySoP/hysop/problem/discrete.py
similarity index 100%
rename from HySoP/hysop/problem/DiscreteProblem.py
rename to HySoP/hysop/problem/discrete.py
diff --git a/HySoP/hysop/problem/problem.py b/HySoP/hysop/problem/problem.py
new file mode 100644
index 000000000..2aa412159
--- /dev/null
+++ b/HySoP/hysop/problem/problem.py
@@ -0,0 +1,90 @@
+# -*- coding: utf-8 -*-
+"""
+@package Problem
+Problem representation
+"""
+
+
+class Problem():
+ """
+ Problem representation.
+ """
+
+ def __init__(self, topology, operators, params=None):
+ """
+ Constructor.
+ Create a transport problem instance.
+
+ """
+ self.domains = []
+ self.variables = []
+ self.operators = operators
+ for op in self.operators:
+ self.addVariable(op.variables)
+ for v in self.variables:
+ self.addDomain(v.domain)
+ self.topology = topology
+ self.solver = None
+
+ def setSolver(self, solver):
+ self.solver = solver(self)
+ ## Discretise domains
+ for d in self.domains:
+ d.discretize(self.topology.resolution)
+ ## Discretise fields
+ for v in self.variables:
+ v.discretize(self.topology)
+ ## Discretise operators
+ for op in self.operators:
+ op.discretize(result_position=self.solver.p_position, result_scalar=self.solver.p_scalar)
+
+ def addVariable(self, cVariable):
+ """
+ Add an continuous variables to the operator.
+ Also add variables' domains to the operator.
+
+ @param cVariable ContinuousVariable.ContinuousVariable : variables to add.
+ """
+ try:
+ for v in cVariable:
+ if self.variables.count(v) == 0:
+ self.variables.append(v)
+ except (TypeError):
+ if self.variables.count(cVariable) == 0:
+ self.variables.append(cVariable)
+
+ def addDomain(self, cDomain):
+ """
+ Add an domain to the operator.
+
+ @param cDomain : domain to add
+ """
+ try:
+ for d in cDomain:
+ if self.domains.count(d) == 0:
+ self.domains.append(d)
+ except (TypeError):
+ if self.domains.count(cDomain) == 0:
+ self.domains.append(cDomain)
+
+ def __str__(self):
+ """ToString method"""
+ s = "Problem based on\n"
+ s += str(self.topology)
+ s += "with following operators : \n"
+ for op in self.operators:
+ s += str(op)
+ s += "working on fields : \n"
+ for f in self.variables:
+ s += str(f)
+ s += "defined on domains : \n"
+ for d in self.domains:
+ s += str(d)
+ s += "\nsolved by : \n" + str(self.solver)
+
+ return s
+
+if __name__ == "__main__":
+ print __doc__
+ print "- Provided class : ContinuousTransportProblem"
+ print ContinuousTransportProblem.__doc__
--
GitLab