From 28ab5d51e9718db002551d6b7d4efe0c45b968ff Mon Sep 17 00:00:00 2001
From: Jonathan Lambrechts <jonathan.lambrechts@uclouvain.be>
Date: Fri, 29 Jan 2010 11:18:04 +0000
Subject: [PATCH] add structured grid function
---
Solver/dgGroupOfElements.cpp | 12 +++--
Solver/function.cpp | 94 ++++++++++++++++++++++++++++++++++--
2 files changed, 99 insertions(+), 7 deletions(-)
diff --git a/Solver/dgGroupOfElements.cpp b/Solver/dgGroupOfElements.cpp
index 8f3d3bafe0..e9a3252442 100644
--- a/Solver/dgGroupOfElements.cpp
+++ b/Solver/dgGroupOfElements.cpp
@@ -337,8 +337,10 @@ dgGroupOfFaces::dgGroupOfFaces (const dgGroupOfElements &elGroup, std::string bo
_groupLeft(elGroup),_groupRight(elGroup)
{
_boundaryTag=boundaryTag;
- if(boundaryTag=="")
- throw;
+ if(boundaryTag==""){
+ Msg::Warning ("empty boundary tag, group of boundary faces not created");
+ return;
+ }
_fsLeft=_groupLeft.getElement(0)->getFunctionSpace(pOrder);
_closuresLeft = _fsLeft->vertexClosure;
_fsRight=NULL;
@@ -357,8 +359,10 @@ dgGroupOfFaces::dgGroupOfFaces (const dgGroupOfElements &elGroup, std::string bo
_groupLeft(elGroup),_groupRight(elGroup)
{
_boundaryTag=boundaryTag;
- if(boundaryTag=="")
- throw;
+ if(boundaryTag==""){
+ Msg::Warning ("empty boundary tag, group of boundary faces not created");
+ return;
+ }
_fsLeft=_groupLeft.getElement(0)->getFunctionSpace(pOrder);
_closuresLeft = _fsLeft->edgeClosure;
_fsRight=NULL;
diff --git a/Solver/function.cpp b/Solver/function.cpp
index 661611d72e..d0f5fbb8e5 100644
--- a/Solver/function.cpp
+++ b/Solver/function.cpp
@@ -1,7 +1,8 @@
+#include <sstream>
+#include <fstream>
#include "function.h"
#include "SPoint3.h"
#include "MElement.h"
-#include <sstream>
// dataCache members
dataCache::dataCache(dataCacheMap *cacheMap) : _valid(false) {
@@ -37,8 +38,10 @@ function *function::get(std::string functionName, bool acceptNull)
if (it==_allFunctions.end()) {
if (acceptNull)
return NULL;
- else
+ else{
+ Msg::Error("unknown function : '%s'\n",functionName.c_str());
throw;
+ }
}
return it->second;
}
@@ -115,6 +118,83 @@ class functionXYZ : public function {
}
};
+class functionStructuredGridFile : public function {
+ public:
+ class data;
+ std::string _coordFunction;
+ int n[3];
+ double d[3],o[3];
+ double get(int i,int j, int k){
+ return v[(i*n[1]+j)*n[2]+k];
+ }
+ double *v;
+ class data : public dataCacheDouble{
+ dataCacheDouble &coord;
+ functionStructuredGridFile *_fun;
+ public:
+ data(functionStructuredGridFile *fun, dataCacheMap *m):
+ dataCacheDouble(*m,m->getNbEvaluationPoints(),1),
+ coord(m->get(fun->_coordFunction,this))
+ {
+ _fun=fun;
+ }
+ void _eval(){
+ for(int pt=0;pt<_value.size1();pt++){
+ double xi[3];
+ int id[3];
+ for(int i=0;i<3;i++){
+ id[i]=(int)((coord(pt,i)-_fun->o[i])/_fun->d[i]);
+ int a=id[i];
+ id[i]=std::max(0,std::min(_fun->n[i]-1,id[i]));
+ xi[i]=(coord(pt,i)-_fun->o[i]-id[i]*_fun->d[i])/_fun->d[i];
+ xi[i]=std::min(1.,std::max(0.,xi[i]));
+ }
+ _value(pt,0) =
+ _fun->get(id[0] ,id[1] ,id[2] )*(1-xi[0])*(1-xi[1])*(1-xi[2])
+ +_fun->get(id[0] ,id[1] ,id[2]+1 )*(1-xi[0])*(1-xi[1])*( xi[2])
+ +_fun->get(id[0] ,id[1]+1 ,id[2] )*(1-xi[0])*( xi[1])*(1-xi[2])
+ +_fun->get(id[0] ,id[1]+1 ,id[2]+1 )*(1-xi[0])*( xi[1])*( xi[2])
+ +_fun->get(id[0]+1 ,id[1] ,id[2] )*( xi[0])*(1-xi[1])*(1-xi[2])
+ +_fun->get(id[0]+1 ,id[1] ,id[2]+1 )*( xi[0])*(1-xi[1])*( xi[2])
+ +_fun->get(id[0]+1 ,id[1]+1 ,id[2] )*( xi[0])*( xi[1])*(1-xi[2])
+ +_fun->get(id[0]+1 ,id[1]+1 ,id[2]+1 )*( xi[0])*( xi[1])*( xi[2]);
+ }
+ }
+ };
+ dataCacheDouble *newDataCache(dataCacheMap* m) {
+ return new data(this,m);
+ }
+ functionStructuredGridFile(const std::string filename, const std::string coordFunction){
+ _coordFunction=coordFunction;
+ std::ifstream input(filename.c_str());
+ if(!input)
+ Msg::Error("cannot open file : %s",filename.c_str());
+ if(filename.substr(filename.size()-4,4)!=".bin") {
+ input>>o[0]>>o[1]>>o[2]>>d[0]>>d[1]>>d[2]>>n[0]>>n[1]>>n[2];
+ int nt = n[0]*n[1]*n[2];
+ v = new double [nt];
+ for (int i=0; i<nt; i++){
+ input>>v[i];
+ }
+ } else {
+ input.read((char *)o, 3 * sizeof(double));
+ input.read((char *)d, 3 * sizeof(double));
+ input.read((char *)n, 3 * sizeof(int));
+ int nt = n[0] * n[1] * n[2];
+ v = new double[nt];
+ input.read((char *)v, nt * sizeof(double));
+ }
+ static int c=0;
+ std::ostringstream oss;
+ oss<<"FunctionStructured"<<c++;
+ _name = oss.str();
+ function::add(_name,this);
+ }
+ ~functionStructuredGridFile(){
+ delete []v;
+ }
+};
+
// constant values copied over each line
class functionConstant::data : public dataCacheDouble {
@@ -180,7 +260,6 @@ public:
};
*/
-
#include "Bindings.h"
void function::registerDefaultFunctions()
@@ -193,11 +272,20 @@ void function::registerBindings(binding *b){
methodBinding *mb;
mb = cb->addMethod("getName",&function::getName);
mb->setDescription("Return the string which identifies this function.");
+
cb = b->addClass<functionConstant>("functionConstant");
cb->setDescription("A constant (scalar or vector) function");
mb = cb->setConstructor<functionConstant,fullMatrix<double>*>();
mb->setArgNames("v",NULL);
mb->setDescription("A new constant function wich values 'v' everywhere. v can be a row-vector.");
cb->setParentClass<function>();
+
+ cb = b->addClass<functionStructuredGridFile>("functionStructuredGridFile");
+ cb->setParentClass<function>();
+ cb->setDescription("A function to interpolate through data given on a structured grid");
+ mb = cb->setConstructor<functionStructuredGridFile,std::string, std::string>();
+ mb->setArgNames("fileName","coordinateFunction",NULL);
+ mb->setDescription("Tri-linearly interpolate through data in file 'fileName' at coordinate given by 'coordinateFunction'.\nThe file format is :\nx0 y0 z0\ndx dy dz\nnx ny nz\nv(0,0,0) v(0,0,1) v(0 0 2) ...");
+
}
--
GitLab