Skip to content
Snippets Groups Projects
Commit 87117cc0 authored by Akash Anand's avatar Akash Anand
Browse files

*** empty log message ***

parent bcc12fc7
No related branches found
No related tags found
No related merge requests found
Hide whitespace changes
Inline Side-by-side
Showing
with 310 additions and 278 deletions
Geo/FEdge.cpp
+ 6
− 62
View file @ 87117cc0
......@@ -4,19 +4,13 @@
#if defined(HAVE_FOURIER_MODEL)
FEdge::FEdge(GModel *model, FM_Edge* edge_, int tag, GVertex *v0, GVertex *v1)
: GEdge(model, tag, v0, v1), edge(edge_), face(0), edgeNum(-1)
: GEdge(model, tag, v0, v1), edge(edge_)
{
//meshAttributes.Method = TRANSFINI;
//meshAttributes.coeffTransfinite = 1.;
//meshAttributes.nbPointsTransfinite = 10;
}
FEdge::FEdge(GModel *model, FM_Face* face_, int edgeNum_, int tag, GVertex *v0,
GVertex *v1)
: GEdge(model, tag, v0, v1), edge(0), face(face_), edgeNum(edgeNum_)
{
}
Range<double> FEdge::parBounds(int i) const
{
return(Range<double>(0.,1.));
......@@ -25,73 +19,23 @@ Range<double> FEdge::parBounds(int i) const
GPoint FEdge::point(double p) const
{
double x, y, z;
if (edge)
edge->F(p,x,y,z);
else {
if (edgeNum == 0) {
double p0, p1, q0, q1;
face->Inverse(v0->x(),v0->y(),v0->z(),p0,q0);
face->Inverse(v1->x(),v1->y(),v1->z(),p1,q1);
face->F(p0+p*(p1-p0),0,x,y,z);
}
else if (edgeNum == 1) {
double p0, p1, q0, q1;
face->Inverse(v0->x(),v0->y(),v0->z(),q0,p0);
face->Inverse(v1->x(),v1->y(),v1->z(),q1,p1);
face->F(1.,p0+p*(p1-p0),x,y,z);
}
else if (edgeNum == 2) {
double p0, p1, q0, q1;
face->Inverse(v0->x(),v0->y(),v0->z(),p0,q0);
face->Inverse(v1->x(),v1->y(),v1->z(),p1,q1);
face->F(p0+p*(p1-p0),1.,x,y,z);
}
else if (edgeNum == 3) {
double p0, p1, q0, q1;
face->Inverse(v0->x(),v0->y(),v0->z(),q0,p0);
face->Inverse(v1->x(),v1->y(),v1->z(),q1,p1);
face->F(0.,p0+p*(p1-p0),x,y,z);
}
else
Msg(INFO,"Invalid edge number.");
}
edge->F(p,x,y,z);
return GPoint(x,y,z);
}
double FEdge::parFromPoint(const SPoint3 &pt) const
{
double p;
if (edge)
edge->Inverse(pt.x(),pt.y(),pt.z(),p);
else {
if ((edgeNum == 0) || (edgeNum == 2)) {
double p0, p1, q0, q1, q;
face->Inverse(v0->x(),v0->y(),v0->z(),p0,q0);
face->Inverse(v1->x(),v1->y(),v1->z(),p1,q1);
face->Inverse(pt.x(),pt.y(),pt.z(),p,q);
p = p0 + p * (p1 - p0);
}
else if ((edgeNum == 1) || (edgeNum == 3)) {
double p0, p1, q0, q1, q;
face->Inverse(v0->x(),v0->y(),v0->z(),q0,p0);
face->Inverse(v1->x(),v1->y(),v1->z(),q1,p1);
face->Inverse(pt.x(),pt.y(),pt.z(),q,p);
p = p0 + p * (p1 - p0);
}
else
Msg(INFO,"Invalid edge number.");
}
edge->Inverse(pt.x(),pt.y(),pt.z(),p);
return p;
}
SVector3 FEdge::firstDer(double par) const
{
double x,y,z;
if (edge)
edge->Dfdt(par,x,y,z);
else {
x = y = z = 0.;
}
edge->Dfdt(par,x,y,z);
return SVector3(x,y,z);
}
......
Geo/FEdge.h
+ 0
− 3
View file @ 87117cc0
......@@ -15,12 +15,9 @@
class FEdge : public GEdge {
protected:
FM_Edge* edge;
FM_Face* face;
int edgeNum;
public:
FEdge(GModel *model, FM_Edge* edge_, int tag, GVertex *v0, GVertex *v1);
FEdge(GModel *model, FM_Face* face_, int edgeNum_, int tag, GVertex *v0,
GVertex *v1);
virtual ~FEdge() {}
double period() const { throw ; }
virtual bool periodic(int dim=0) const { return false; }
......
Geo/FFace.cpp
+ 1
− 1
View file @ 87117cc0
......@@ -58,7 +58,7 @@ int FFace::containsParam(const SPoint2 &pt) const
if(pt[1] < 0. - tol || pt[1] > 1. + tol) return 0;
return 1;
}
SVector3 FFace::normal(const SPoint2 &param) const
{
double x,y,z;
......
Mesh/meshGEdge.cpp
+ 2
− 1
View file @ 87117cc0
// $Id: meshGEdge.cpp,v 1.38 2007-05-02 07:59:27 geuzaine Exp $
// $Id: meshGEdge.cpp,v 1.39 2007-05-12 22:45:14 anand Exp $
//
// Copyright (C) 1997-2007 C. Geuzaine, J.-F. Remacle
//
......@@ -227,6 +227,7 @@ void meshGEdge::operator() (GEdge *ge)
a = Integration(ge, t_begin, t_end, F_Lc, Points, 1.e-8);
N = std::max(ge->minimumMeshSegments() + 1, (int)(a + 1.));
}
const double b = a / (double)(N - 1);
// if the curve is periodic and if the begin vertex is identical to the end vertex
......
contrib/FourierModel/ContinuationPatch.cpp
+ 2
− 2
View file @ 87117cc0
......@@ -166,8 +166,8 @@ std::complex<double> ContinuationPatch::
double epsilon = 1e-12;
if(u < 0. - epsilon || u > 1. + epsilon){
Msg::Error("Trying to interpolate outside interval: (u,v)=(%.16g,%.16g) "
"not in [%g,%g]x[%g,%g]", u, v, 0., 1., 0., 1.);
}
"not in [%g,%g]x[%g,%g]", u, v, 0., 1., 0., 1.);
}
// Interpolate to find value at (u,v)
for(int j = 0; j < _uModes; j++){
......
contrib/FourierModel/Curve.cpp
+ 271
− 192
View file @ 87117cc0
......@@ -2,228 +2,307 @@
#include "Message.h"
#include "Curve.h"
Curve::Curve(IntersectionInfo* II, std::vector<Patch*> patches) : _II(II)
Curve::Curve(IntersectionInfo* II, std::vector<Patch*> patches)
: _II(II)
{
_h = 1.e-6;
_tol = 1.e-12;
for (int i=0;i<patches.size();i++) {
if (patches[i]->GetTag() == _II->intersectingPatches[0].patchTag) {
_patch0 = patches[i];
break;
if (_II->intersectingPatches[0].patchTag < 0)
_patch0 = 0;
else {
for (int i=0;i<patches.size();i++) {
if (patches[i]->GetTag() == _II->intersectingPatches[0].patchTag) {
_patch0 = patches[i];
break;
}
}
}
for (int i=0;i<patches.size();i++)
if (patches[i]->GetTag() == _II->intersectingPatches[1].patchTag) {
_patch1 = patches[i];
break;
}
/*
int n1=64;
int n2=64;
double h1 = 1000./(n1-1);
double h2 = 1000./(n2-1);
for (int j=0;j<n1;j++) {
for (int k=0;k<n2;k++) {
double u = j*h1;
double v = k*h2;
double x, y, z;
_patch1->F(u,v,x,y,z);
printf("x(%d,%d) = %g; y(%d,%d) = %g; z(%d,%d) = %g;\n",
j+1,k+1,x,j+1,k+1,y,j+1,k+1,z);
if (_II->intersectingPatches[1].patchTag)
_patch1 = 0;
else {
for (int i=0;i<patches.size();i++) {
if (patches[i]->GetTag() == _II->intersectingPatches[1].patchTag) {
_patch1 = patches[i];
break;
}
}
}
*/
double u,v,x,y,z;
x = _II->SP[0]; y = _II->SP[1]; z = _II->SP[2];
_patch0->Inverse(x,y,z,u,v);
_SP0[0] = u;
_SP0[1] = v;
_patch1->Inverse(x,y,z,u,v);
_SP1[0] = u;
_SP1[1] = v;
//printf("SPx = %g; SPy = %g; SPz = %g;\n",x,y,z);
x = _II->EP[0]; y = _II->EP[1]; z = _II->EP[2];
_patch0->Inverse(x,y,z,u,v);
_EP0[0] = u;
_EP0[1] = v;
_patch1->Inverse(x,y,z,u,v);
_EP1[0] = u;
_EP1[1] = v;
//printf("EPx = %g; EPy = %g; EPz = %g;\n",x,y,z);
if ((_patch0) && (_patch1)) {
double u,v,x,y,z;
x = _II->SP[0]; y = _II->SP[1]; z = _II->SP[2];
_patch0->Inverse(x,y,z,u,v);
_SP0[0] = u;
_SP0[1] = v;
_patch1->Inverse(x,y,z,u,v);
_SP1[0] = u;
_SP1[1] = v;
//printf("SPx = %g; SPy = %g; SPz = %g;\n",x,y,z);
x = _II->EP[0]; y = _II->EP[1]; z = _II->EP[2];
_patch0->Inverse(x,y,z,u,v);
_EP0[0] = u;
_EP0[1] = v;
_patch1->Inverse(x,y,z,u,v);
_EP1[0] = u;
_EP1[1] = v;
//printf("EPx = %g; EPy = %g; EPz = %g;\n",x,y,z);
}
else if (_patch0) {
double u,v,x,y,z;
x = _II->SP[0]; y = _II->SP[1]; z = _II->SP[2];
_patch0->Inverse(x,y,z,u,v);
_SP0[0] = u;
_SP0[1] = v;
x = _II->EP[0]; y = _II->EP[1]; z = _II->EP[2];
_patch0->Inverse(x,y,z,u,v);
_EP0[0] = u;
_EP0[1] = v;
}
else if (_patch1) {
double u,v,x,y,z;
x = _II->SP[0]; y = _II->SP[1]; z = _II->SP[2];
_patch1->Inverse(x,y,z,u,v);
_SP1[0] = u;
_SP1[1] = v;
x = _II->EP[0]; y = _II->EP[1]; z = _II->EP[2];
_patch1->Inverse(x,y,z,u,v);
_EP1[0] = u;
_EP1[1] = v;
}
else
Msg::Error("No patches to find the intersection curve.");
}
void Curve::F(double t, double &x, double &y, double &z)
{
double u0, v0, u1, v1;
//Msg::Info("patch0 : %d",_patch0->GetTag());
//Msg::Info("patch1 : %d",_patch1->GetTag());
//Msg::Info("%g,%g,%g",_II->SP[0],_II->SP[1],_II->SP[2]);
//Msg::Info("%g,%g,%g",_II->EP[0],_II->EP[1],_II->EP[2]);
//Msg::Info("%g,%g,%g,%g",_SP0[0],_SP0[1],_EP0[0],_EP0[1]);
//Msg::Info("%g,%g,%g,%g",_SP1[0],_SP1[1],_EP1[0],_EP1[1]);
if ((std::abs(_SP0[0]-_EP0[0])<1.e-12) &&
(_patch0->_PI->periodic[0] == 1)) {
u0 = _SP0[0] + t * (1. + _EP0[0] - _SP0[0]);
if (u0 > 1.)
if ((_patch0) && (_patch1)) {
double u0, v0, u1, v1;
//printf("patch0 : %d\n",_patch0->GetTag());
//printf("patch1 : %d\n",_patch1->GetTag());
//printf("%g,%g,%g\n",_II->SP[0],_II->SP[1],_II->SP[2]);
//printf("%g,%g,%g\n",_II->EP[0],_II->EP[1],_II->EP[2]);
//printf("%g,%g,%g,%g\n",_SP0[0],_SP0[1],_EP0[0],_EP0[1]);
//printf("%g,%g,%g,%g\n",_SP1[0],_SP1[1],_EP1[0],_EP1[1]);
if (_patch0->_PI->periodic[0] == 1) {
if (std::abs(_SP0[0]-_EP0[0])<1.e-12)
u0 = _SP0[0] + t * (1. + _EP0[0] - _SP0[0]);
else
u0 = _SP0[0] + t * (_EP0[0] - _SP0[0]);
u0 -= floor(u0);
}
else
u0 = _SP0[0] + t * (_EP0[0] - _SP0[0]);
if ((std::abs(_SP0[1]-_EP0[1])<1.e-12) &&
(_patch0->_PI->periodic[1] == 1)) {
v0 = _SP0[1] + t * (1. + _EP0[1] - _SP0[1]);
if (u0 > 1.)
}
else
u0 = _SP0[0] + t * (_EP0[0] - _SP0[0]);
if (_patch0->_PI->periodic[1] == 1) {
if (std::abs(_SP0[1]-_EP0[1])<1.e-12)
v0 = _SP0[1] + t * (1. + _EP0[1] - _SP0[1]);
else
v0 = _SP0[1] + t * (_EP0[1] - _SP0[1]);
v0 -= floor(v0);
}
else
v0 = _SP0[1] + t * (_EP0[1] - _SP0[1]);
if ((std::abs(_SP1[0]-_EP1[0])<1.e-12) &&
(_patch1->_PI->periodic[0] == 1)) {
u1 = _SP1[0] + t * (1. + _EP1[0] - _SP1[0]);
if (u1 > 1.)
}
else
v0 = _SP0[1] + t * (_EP0[1] - _SP0[1]);
if (_patch1->_PI->periodic[0] == 1) {
if (std::abs(_SP1[0]-_EP1[0])<1.e-12)
u1 = _SP1[0] + t * (1. + _EP1[0] - _SP1[0]);
else
u1 = _SP1[0] + t * (_EP1[0] - _SP1[0]);
u1 -= floor(u1);
}
else {
u1 = _SP1[0] + t * (_EP1[0] - _SP1[0]);
}
if ((std::abs(_SP1[1]-_EP1[1])<1.e-12) &&
(_patch1->_PI->periodic[1] == 1)) {
v1 = _SP1[1] + t * (1. + _EP1[1] - _SP1[1]);
if (u1 > 1.)
}
else
u1 = _SP1[0] + t * (_EP1[0] - _SP1[0]);
if (_patch1->_PI->periodic[1] == 1) {
if (std::abs(_SP1[1]-_EP1[1])<1.e-12)
v1 = _SP1[1] + t * (1. + _EP1[1] - _SP1[1]);
else
v1 = _SP1[1] + t * (_EP1[1] - _SP1[1]);
v1 -= floor(v1);
}
else
v1 = _SP1[1] + t * (_EP1[1] - _SP1[1]);
//Msg::Info("%g,%g,%g,%g",u0,v0,u1,v1);
double x0, y0, z0;
_patch0->F(u0,v0,x0,y0,z0);
//Msg::Info("%g,%g,%g",x0,y0,z0);
double x1, y1, z1;
_patch1->F(u1,v1,x1,y1,z1);
//Msg::Info("%g,%g,%g",x1,y1,z1);
double r,u,v,rPlus,uPlus,vPlus;
r = u0;
u = u1;
v = v1;
double F[3];
F[0] = x1 - x0;
F[1] = y1 - y0;
F[2] = z1 - z0;
while ((std::abs(F[0])>_tol) || (std::abs(F[1])>_tol) ||
(std::abs(F[2])>_tol)) {
rPlus = r + _h;
uPlus = u + _h;
vPlus = v + _h;
if (_patch0->_PI->periodic[0] == 1)
rPlus -= std::floor(rPlus);
if (_patch1->_PI->periodic[0] == 1)
uPlus -= std::floor(uPlus);
if (_patch1->_PI->periodic[1] == 1)
vPlus -= floor(vPlus);
//Msg::Info(" F : %g,%g,%g",F[0],F[1],F[2]);
//Msg::Info(" R : %g,%g,%g",x0,y0,z0);
double x0rPlus, y0rPlus, z0rPlus;
_patch0->F(rPlus,v0,x0rPlus,y0rPlus,z0rPlus);
double x1uPlus, y1uPlus, z1uPlus;
_patch1->F(uPlus,v,x1uPlus,y1uPlus,z1uPlus);
double x1vPlus, y1vPlus, z1vPlus;
_patch1->F(u,vPlus,x1vPlus,y1vPlus,z1vPlus);
}
else
v1 = _SP1[1] + t * (_EP1[1] - _SP1[1]);
double Df[3][3];
Df[0][0] = - (x0rPlus - x0) / _h;
Df[0][1] = (x1uPlus - x1) / _h;
Df[0][2] = (x1vPlus - x1) / _h;
Df[1][0] = - (y0rPlus - y0) / _h;
Df[1][1] = (y1uPlus - y1) / _h;
Df[1][2] = (y1vPlus - y1) / _h;
Df[2][0] = - (z0rPlus - z0) / _h;
Df[2][1] = (z1uPlus - z1) / _h;
Df[2][2] = (z1vPlus - z1) / _h;
//Msg::Info("D1: %g,%g,%g",Df[0][0],Df[0][1],Df[0][2]);
//Msg::Info("D2: %g,%g,%g",Df[1][0],Df[1][1],Df[1][2]);
//Msg::Info("D3: %g,%g,%g",Df[2][0],Df[2][1],Df[2][2]);
//printf("%d,%d,%g,%g: %g,%g,%g,%g\n",_patch1->_PI->periodic[0],
// _patch1->_PI->periodic[1],std::abs(_SP1[0]-_EP1[0]),
// std::abs(_SP1[1]-_EP1[1]),u0,v0,u1,v1);
double det =
Df[0][0] * (Df[1][1] * Df[2][2] - Df[1][2] * Df[2][1]) +
Df[0][1] * (Df[1][2] * Df[2][0] - Df[1][0] * Df[2][2]) +
Df[0][2] * (Df[1][0] * Df[2][1] - Df[1][1] * Df[2][0]);
double x0, y0, z0;
_patch0->F(u0,v0,x0,y0,z0);
//Msg::Info("det: %g",det);
double update[3];
update[0] =
(Df[1][1] * Df[2][2] - Df[1][2] * Df[2][1]) * F[0] +
(Df[0][2] * Df[2][1] - Df[0][1] * Df[2][2]) * F[1] +
(Df[0][1] * Df[1][2] - Df[0][2] * Df[1][1]) * F[2];
update[1] =
(Df[1][2] * Df[2][0] - Df[1][0] * Df[2][2]) * F[0] +
(Df[0][0] * Df[2][2] - Df[0][2] * Df[2][0]) * F[1] +
(Df[0][2] * Df[1][0] - Df[0][0] * Df[1][2]) * F[2];
update[2] =
(Df[1][0] * Df[2][1] - Df[1][1] * Df[2][0]) * F[0] +
(Df[0][1] * Df[2][0] - Df[0][0] * Df[2][1]) * F[1] +
(Df[0][0] * Df[1][1] - Df[0][1] * Df[1][0]) * F[2];
//Msg::Info("U: %g,%g,%g",update[0],update[1],update[2]);
r -= update[0] / det;
u -= update[1] / det;
v -= update[2] / det;
//Msg::Info("U/det: %g,%g,%g",update[0]/det,update[1]/det,update[2]/det);
//Msg::Info("uv: %g,%g,%g",r,u,v);
if (_patch0->_PI->periodic[0] == 1)
r -= std::floor(r);
if (_patch1->_PI->periodic[0] == 1)
u -= std::floor(u);
if (_patch1->_PI->periodic[1] == 1)
v -= floor(v);
//Msg::Info("UV : %g,%g,%g",r,u,v);
//rPlus = r + _h;
//uPlus = u + _h;
//vPlus = v + _h;
//Msg::Info("%g,%g,%g",x0,y0,z0);
_patch0->F(r,v0,x0,y0,z0);
_patch1->F(u,v,x1,y1,z1);
double x1, y1, z1;
_patch1->F(u1,v1,x1,y1,z1);
//Msg::Info("%g,%g,%g",x1,y1,z1);
double r,u,v,rPlus,uPlus,vPlus;
r = u0;
u = u1;
v = v1;
double F[3];
F[0] = x1 - x0;
F[1] = y1 - y0;
F[2] = z1 - z0;
while ((std::abs(F[0])>_tol) || (std::abs(F[1])>_tol) ||
(std::abs(F[2])>_tol)) {
rPlus = r + _h;
uPlus = u + _h;
vPlus = v + _h;
if (_patch0->_PI->periodic[0] == 1)
rPlus -= std::floor(rPlus);
if (_patch1->_PI->periodic[0] == 1)
uPlus -= std::floor(uPlus);
if (_patch1->_PI->periodic[1] == 1)
vPlus -= floor(vPlus);
//Msg::Info(" F : %g,%g,%g",F[0],F[1],F[2]);
//Msg::Info(" R : %g,%g,%g",x0,y0,z0);
double x0rPlus, y0rPlus, z0rPlus;
_patch0->F(rPlus,v0,x0rPlus,y0rPlus,z0rPlus);
double x1uPlus, y1uPlus, z1uPlus;
_patch1->F(uPlus,v,x1uPlus,y1uPlus,z1uPlus);
double x1vPlus, y1vPlus, z1vPlus;
_patch1->F(u,vPlus,x1vPlus,y1vPlus,z1vPlus);
double Df[3][3];
Df[0][0] = - (x0rPlus - x0) / _h;
Df[0][1] = (x1uPlus - x1) / _h;
Df[0][2] = (x1vPlus - x1) / _h;
Df[1][0] = - (y0rPlus - y0) / _h;
Df[1][1] = (y1uPlus - y1) / _h;
Df[1][2] = (y1vPlus - y1) / _h;
Df[2][0] = - (z0rPlus - z0) / _h;
Df[2][1] = (z1uPlus - z1) / _h;
Df[2][2] = (z1vPlus - z1) / _h;
//Msg::Info("D1: %g,%g,%g",Df[0][0],Df[0][1],Df[0][2]);
//Msg::Info("D2: %g,%g,%g",Df[1][0],Df[1][1],Df[1][2]);
//Msg::Info("D3: %g,%g,%g",Df[2][0],Df[2][1],Df[2][2]);
double det =
Df[0][0] * (Df[1][1] * Df[2][2] - Df[1][2] * Df[2][1]) +
Df[0][1] * (Df[1][2] * Df[2][0] - Df[1][0] * Df[2][2]) +
Df[0][2] * (Df[1][0] * Df[2][1] - Df[1][1] * Df[2][0]);
//Msg::Info("det: %g",det);
double update[3];
update[0] =
(Df[1][1] * Df[2][2] - Df[1][2] * Df[2][1]) * F[0] +
(Df[0][2] * Df[2][1] - Df[0][1] * Df[2][2]) * F[1] +
(Df[0][1] * Df[1][2] - Df[0][2] * Df[1][1]) * F[2];
update[1] =
(Df[1][2] * Df[2][0] - Df[1][0] * Df[2][2]) * F[0] +
(Df[0][0] * Df[2][2] - Df[0][2] * Df[2][0]) * F[1] +
(Df[0][2] * Df[1][0] - Df[0][0] * Df[1][2]) * F[2];
update[2] =
(Df[1][0] * Df[2][1] - Df[1][1] * Df[2][0]) * F[0] +
(Df[0][1] * Df[2][0] - Df[0][0] * Df[2][1]) * F[1] +
(Df[0][0] * Df[1][1] - Df[0][1] * Df[1][0]) * F[2];
//Msg::Info("U: %g,%g,%g",update[0],update[1],update[2]);
r -= update[0] / det;
u -= update[1] / det;
v -= update[2] / det;
//Msg::Info("U/det: %g,%g,%g",update[0]/det,update[1]/det,update[2]/det);
//Msg::Info("uv: %g,%g,%g",r,u,v);
if (_patch0->_PI->periodic[0] == 1)
r -= std::floor(r);
if (_patch1->_PI->periodic[0] == 1)
u -= std::floor(u);
if (_patch1->_PI->periodic[1] == 1)
v -= floor(v);
//Msg::Info("UV : %g,%g,%g",r,u,v);
//rPlus = r + _h;
//uPlus = u + _h;
//vPlus = v + _h;
_patch0->F(r,v0,x0,y0,z0);
_patch1->F(u,v,x1,y1,z1);
F[0] = x1 - x0;
F[1] = y1 - y0;
F[2] = z1 - z0;
}
x = x0; y = y0; z = z0;
}
else if (_patch0) {
if (_II->edgeNumber == 0)
_patch0->F(t,0.,x,y,z);
else if (_II->edgeNumber == 1)
_patch0->F(1.,t,x,y,z);
else if (_II->edgeNumber == 2)
_patch0->F(1.-t,1.,x,y,z);
else if (_II->edgeNumber == 3)
_patch0->F(0.,1.-t,x,y,z);
else
Msg::Error("Unknown edge number.");
}
x = x0; y = y0; z = z0;
else if (_patch1) {
if (_II->edgeNumber == 0)
_patch1->F(t,0.,x,y,z);
else if (_II->edgeNumber == 1)
_patch1->F(1.,t,x,y,z);
else if (_II->edgeNumber == 2)
_patch1->F(1.-t,1.,x,y,z);
else if (_II->edgeNumber == 3)
_patch1->F(0.,1.-t,x,y,z);
}
else
Msg::Error("No patches to find the intersection curve.");
}
bool Curve::Inverse(double x,double y,double z,double &t)
{
double u0, v0;
_patch0->Inverse(x,y,z,u0,v0);
if ((std::abs(_SP0[1]-_EP0[1])<1.e-12) &&
(_patch0->_PI->periodic[1] == 1)) {
t = (v0 - _SP0[1] > 0 ? v0 - _SP0[1] : 1. + v0 - _SP0[1]);
if ((_patch0) && (_patch1)) {
double u0, v0;
_patch0->Inverse(x,y,z,u0,v0);
if ((std::abs(_SP0[1]-_EP0[1])<1.e-12) &&
(_patch0->_PI->periodic[1] == 1)) {
t = (v0 - _SP0[1] > 0 ? v0 - _SP0[1] : 1. + v0 - _SP0[1]);
}
else
t = (v0 - _SP0[1]) / (_EP0[1] - _SP0[1]);
}
else if (_patch0) {
double u,v;
_patch0->Inverse(x,y,z,u,v);
if (_II->edgeNumber == 0)
t = u;
else if (_II->edgeNumber == 1)
t = v;
else if (_II->edgeNumber == 2)
t = 1. - u;
else if (_II->edgeNumber == 3)
t = 1. - v;
}
else if (_patch1) {
double u,v;
_patch1->Inverse(x,y,z,u,v);
if (_II->edgeNumber == 0)
t = u;
else if (_II->edgeNumber == 1)
t = v;
else if (_II->edgeNumber == 2)
t = 1. - u;
else if (_II->edgeNumber == 3)
t = 1. - v;
}
else
t = (v0 - _SP0[1]) / (_EP0[1] - _SP0[1]);
Msg::Error("No patches to find the intersection curve.");
if ((t < 0.) || (t > 1.))
return false;
......
contrib/FourierModel/Curve.h
+ 1
− 0
View file @ 87117cc0
......@@ -7,6 +7,7 @@
// The base class for the patches
class Curve {
private:
int _edgeNum;
double _h, _tol;
protected:
// Patches
......
contrib/FourierModel/FM_Edge.cpp
+ 12
− 12
View file @ 87117cc0
......@@ -76,19 +76,19 @@ void FM_Edge::Dfdt(double t, double &x, double &y, double &z)
}
else
tRescaled = tStart + t * (tEnd - tStart);
if (t+0.5*h > 1.) {
if (tRescaled+0.5*h > tEnd) {
_curve->F(tRescaled, xPlus, yPlus, zPlus);
double tMinus = tStart + (t - h) * (tEnd - tStart);
double tMinus = tRescaled - h;
_curve->F(tMinus, xMinus, yMinus, zMinus);
}
else if (t-0.5*h < 0.) {
else if (tRescaled-0.5*h < tStart) {
_curve->F(tRescaled, xMinus, yMinus, zMinus);
double tPlus = tStart + (t + h) * (tEnd - tStart);
double tPlus = tRescaled + h;
_curve->F(tPlus, xPlus, yPlus, zPlus);
}
else {
double tPlus = tStart + (t + 0.5 * h) * (tEnd - tStart);
double tMinus = tStart + (t - 0.5 * h) * (tEnd - tStart);
double tPlus = tRescaled + 0.5 * h;
double tMinus = tRescaled - 0.5 * h;
_curve->F(tPlus, xPlus, yPlus, zPlus);
_curve->F(tMinus, xMinus, yMinus, zMinus);
}
......@@ -121,19 +121,19 @@ void FM_Edge::Dfdfdtdt(double t, double &x, double &y, double &z)
}
else
tRescaled = tStart + t * (tEnd - tStart);
if (t+0.5*h > 1.) {
if (tRescaled+0.5*h > tEnd) {
Dfdt(tRescaled, xPlus, yPlus, zPlus);
double tMinus = tStart + (t - h) * (tEnd - tStart);
double tMinus = tRescaled - h;
Dfdt(tMinus, xMinus, yMinus, zMinus);
}
else if (t-0.5*h < 0.) {
else if (tRescaled-0.5*h < tStart) {
Dfdt(tRescaled, xMinus, yMinus, zMinus);
double tPlus = tStart + (t + h) * (tEnd - tStart);
double tPlus = tRescaled + h;
Dfdt(tPlus, xPlus, yPlus, zPlus);
}
else {
double tPlus = tStart + (t + 0.5 * h) * (tEnd - tStart);
double tMinus = tStart + (t - 0.5 * h) * (tEnd - tStart);
double tPlus = tRescaled + 0.5 * h;
double tMinus = tRescaled - 0.5 * h;
Dfdt(tPlus, xPlus, yPlus, zPlus);
Dfdt(tMinus, xMinus, yMinus, zMinus);
}
......
contrib/FourierModel/FM_Info.h
+ 2
− 1
View file @ 87117cc0
......@@ -29,8 +29,9 @@ class IntersectionInfo {
struct {
int patchTag;
} intersectingPatches[2];
int edgeNumber;
IntersectionInfo() {}
IntersectionInfo() : tag(-1), edgeNumber(-1) {}
virtual ~IntersectionInfo() {}
};
......
contrib/FourierModel/FM_Reader.cpp
+ 13
− 4
View file @ 87117cc0
......@@ -42,6 +42,7 @@ FM_Reader::FM_Reader(const char* fn)
}
}
}
for (int i=0;i<_nIntersections;i++) {
_intersectionList.push_back(new IntersectionInfo);
InputFile >> _intersectionList[i]->tag;
......@@ -51,10 +52,13 @@ FM_Reader::FM_Reader(const char* fn)
>> _intersectionList[i]->EP[2];
InputFile >> _intersectionList[i]->intersectingPatches[0].patchTag;
InputFile >> _intersectionList[i]->intersectingPatches[1].patchTag;
if ((_intersectionList[i]->intersectingPatches[0].patchTag < 0) ||
(_intersectionList[i]->intersectingPatches[1].patchTag < 0))
InputFile >> _intersectionList[i]->edgeNumber;
}
_patch.resize(_nPatches);
_intersection.resize(_nIntersections);
_intersection.resize(_nIntersections,0);
for (int i=0;i<_nPatches;i++) {
PatchInfo* PI = _patchList[i];
......@@ -109,7 +113,12 @@ Patch* FM_Reader::GetPatch(int tag)
Curve* FM_Reader::GetIntersection(int tag)
{
for (int i=0;i<_intersection.size();i++)
if (_intersection[i]->GetTag() == tag)
return _intersection[i];
Curve* curve = 0;
for (int i=0;i<_intersection.size();i++) {
if (_intersection[i]->GetTag() == tag) {
curve = _intersection[i];
break;
}
}
return curve;
}
0% or .
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment