Commit 1160d11e authored by Phillip Lino Rall's avatar Phillip Lino Rall
Browse files

face transform now working FINE ! Implementation is not perfect though, but...

face transform now working FINE ! Implementation is not perfect though, but since this is only called once, shouldn't make a big difference in performoance.
parent 09102163
......@@ -57,6 +57,7 @@ bool compareVectors(D3vector v, D3vector w)
}
else
{
// cout << "fine " << endl;
return false;
}
}
......@@ -77,7 +78,7 @@ Blockgrid::Blockgrid ( Unstructured_grid *ug_ ) {
number_points[i] = 10;
}
variable_set = false;
}
}
Blockgrid::Blockgrid ( Unstructured_grid *ug_, int N ) {
id_of_grid = id_count_grid; ++id_count_grid;
......@@ -354,6 +355,7 @@ inline std::vector<int> find_p_quad ( dir3D quad, Hexahedron_el *hex )
// WDed, EDed, WTed, ETed,
// SWed, SEed, NWed, NEed };
std::vector<int> ijk = {1 , 1 , 1};
if (quad == Wdir3D)
{
hex->orientation[4];
......@@ -511,6 +513,12 @@ D3vector Blockgrid::Give_coord_hexahedron ( int id_hex,
std::vector<int> p_quad_Ddir3D = find_p_quad(Ddir3D,hex);
std::vector<int> p_quad_Tdir3D = find_p_quad(Tdir3D,hex);
// cout << "obacht, sensitive changes " << endl;
// p_quad_Edir3D[0] = 1;
// p_quad_Edir3D[1] = 1;
int idW = hex->Give_id_quadrangle(Wdir3D);
int idE = hex->Give_id_quadrangle(Edir3D);
int idS = hex->Give_id_quadrangle(Sdir3D);
......@@ -534,7 +542,8 @@ D3vector Blockgrid::Give_coord_hexahedron ( int id_hex,
D3vector QPsi[6], PresQuad;
if (transformFromQuadrangle)
//if (transformFromQuadrangle)
if(ug->Give_transform_From_Quadrangle())
{
//geht bisher nur fuer 1 hexaeder, allgemien so richtig?
// idW = 0 ; idE = 1 ; idS = 2 ; idN = 3 ; idD = 4 ; idT = 5;
......@@ -654,6 +663,8 @@ D3vector Blockgrid::Give_coord_hexahedron ( int id_hex,
//compareVectors(PNW,Give_coord_quadrangle(idS,i,Give_Nz_hexahedron ( id_hex ))); //yes
// compareVectors(PNW,Give_coord_quadrangle(idS,quadSaddI+i*p_quad_Sdir3D[0] ,Give_Nz_hexahedron ( id_hex )-quadSaddK)); //modified
//compareVectors(PNW,Give_coord_quadrangle(idT,i,0)); //no
//PNW = Give_coord_quadrangle(idS,i,Give_Nz_hexahedron ( id_hex )); //ok
PNW = Give_coord_quadrangle(idS,quadSaddI+i*p_quad_Sdir3D[0] ,nEdgeZ-quadSaddK);
......@@ -989,75 +1000,174 @@ D3vector Blockgrid::Give_coord_edge ( int id, int i ) const
{
//Phillip hinzu
//Funktioniert jetzt, aber super ineffizient. quadFix sollte gespeichert werden, da entlang einer kante konstant.
//--> optimierungsbedarf, aber erstmal OK.
//ACHTUNG : an kanten stimmt es noch nicht ganz. ähnlich wie bei give_coord_hexahedra, manchmal stimmt die richtung nicht ( n-i anstatt i) oder n anstatt 0.
//test überlegen, der das hinbekommt? vllt einfacher:
Quadrangle_el* quad;
Hexahedron_el* hex;
D3vector Psi;
//if ()
for (int quadIter = 0; quadIter < ug->Give_number_quadrangles();quadIter++)
if (ug->Give_transform_From_Quadrangle())
{
int idHex;
dir3D dirForQuad;
int quadId;
if (ug->Give_edge_to_quad_id(id) == -1)
{
// cout << "slow way, once each, id " << id << endl;
//findet hexahedron entsprechend zur kanten ID
for (int iterHex = 0 ; iterHex < ug->Give_number_hexahedra() ; iterHex++)
{
idHex = iterHex;
hex = ug->Give_hexahedron(iterHex);
for (int iter = 0 ; iter < 12 ; iter++)
{
int idCorner = hex->Give_id_edge(Edges_cell(iter));
if (idCorner == id)
{
iterHex = ug->Give_number_hexahedra() ;
iter = 12;
}
}
}
//findet quad ID sowie ausrichtung entsprechend zur kanten ID
for ( int iter = 0 ; iter < 6 ; iter++)
{
quad = ug->Give_quadrangle(hex->Give_id_quadrangle((dir3D)iter));
dirForQuad = (dir3D)iter ;
for (int iterEdge = 0 ; iterEdge < 4 ; iterEdge++)
{
quadId = quad->Give_id_edge((dir2D)iterEdge);
if (id == quadId)
{
quadId = hex->Give_id_quadrangle((dir3D)iter);
iterEdge = 4;
iter = 6;
}
}
}
ug->Set_edge_to_quad_id(id,quadId);
ug->Set_edge_to_quad_dir(id,(int)dirForQuad);
}
else
{
// dirForQuad =
// cout << "fast way, after init only, id " << id << endl;
quadId = ug->Give_edge_to_quad_id(id);
dirForQuad = (dir3D)ug->Give_edge_to_quad_dir(id);
quad = ug->Give_quadrangle(quadId);
}
int Nx = Give_Nx_quadrangle(dirForQuad);
int Ny = Give_Ny_quadrangle(dirForQuad);
// int quadAddX = (quadFix[0] == 1) ? 0 : Nx;
// int quadAddY = (quadFix[1] == 1) ? 0 : Ny;
//cout << "quad ids: " << quad->Give_id_edge(Wdir2D) << " " << quad->Give_id_edge(Edir2D) << " " << quad->Give_id_edge(Ndir2D) << " " << quad->Give_id_edge(Sdir2D) << endl;
if (id == quad->Give_id_edge(Wdir2D)) //Wdir2D
{
// cout << "Wdir2D" << endl;
//Psi = Give_coord_quadrangle(quadId,0,i);
Psi = Give_coord_quadrangle(quadId,0,i);
// D3vector Psi2 = Give_coord_quadrangle(quadId,quadAddX,quadAddY+i*quadFix[1]);
// compareVectors(Psi,Psi2);
// Psi = Give_coord_quadrangle(quadId,quadAddX,quadAddY+i*quadFix[1]);
}
if (id == quad->Give_id_edge(Edir2D)) //Edir2D
{
// cout << "Edir2D" << endl;
// Psi = Give_coord_quadrangle(quadId,Give_Nx_hexahedron ( 0 ),i);
Psi = Give_coord_quadrangle(quadId,Nx,i);
// D3vector Psi2 = Give_coord_quadrangle(quadId,Nx-quadAddX,quadAddY+i*quadFix[1]);
// compareVectors(Psi,Psi2);
// Psi = Give_coord_quadrangle(quadId,Nx-quadAddX,quadAddY+i*quadFix[1]);
}
if (id == quad->Give_id_edge(Sdir2D)) //Sdir2D
{
// cout << "Sdir2D" << endl;
// Psi = Give_coord_quadrangle(quadId,i,0);
Psi = Give_coord_quadrangle(quadId,i,0);
// D3vector Psi2 = Give_coord_quadrangle(quadId,quadAddX+i*quadFix[0],quadAddY);
// compareVectors(Psi,Psi2);
// Psi = Give_coord_quadrangle(quadId,quadAddX+i*quadFix[0],quadAddY);
}
if (id == quad->Give_id_edge(Ndir2D)) //Ndir2D
{
// cout << "Ndir2D" << endl;
// Psi = Give_coord_quadrangle(quadId,i,Give_Nx_hexahedron ( 0 ));
Psi = Give_coord_quadrangle(quadId,i,Ny);
// D3vector Psi2 = Give_coord_quadrangle(quadId,quadAddX+i*quadFix[0],Ny - quadAddY);
// compareVectors(Psi,Psi2);
// Psi = Give_coord_quadrangle(quadId,quadAddX+i*quadFix[0],Ny - quadAddY);
}
// if (fabs(sqrt(Psi.x*Psi.x+Psi.y*Psi.y +Psi.z*Psi.z ) - 10) < 0.01 && Psi.z > -8.660254037844387)
// {
// cout << "quadId " << quadId << endl;
// cout << "dirForQuad " << dirForQuad << endl;
// cout << "idHex " << idHex << endl;
// cout << "Psi.x " << Psi.x << endl;
// cout << "Psi.y " << Psi.y << endl;
// cout << "xx+yy" << Psi.x*Psi.x+Psi.y*Psi.y << endl;
// cout << "Psi.z " << Psi.y << endl;
// cout << endl;
// }
return Psi;
//new version !? maybe faster
for (int quadIter = 0; quadIter < ug->Give_number_quadrangles();quadIter++)
{
quad = ug->Give_quadrangle(quadIter);
if (id == quad->Give_id_edge(Wdir2D))
{
// cout << "Quadrangle id: " << quadIter << " ";
// cout << "correct edge found: Wdir2D " << endl;
// Give_coord_quadrangle(quadIter,i,0).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,0,i).Print();cout<<endl; // korrekt
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 )-i,0).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,0,Give_Nx_hexahedron ( 0 )-i).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,i,Give_Nx_hexahedron ( 0 )).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 ),i).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 )-i,Give_Nx_hexahedron ( 0 )).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 ),Give_Nx_hexahedron ( 0 )-i).Print();cout<<endl;
Psi = Give_coord_quadrangle(quadIter,0,i);
//D3vector Psi2 = Give_coord_quadrangle(quadIter,quadAddX,quadAddY+i*quadFix[1]);
//compareVectors(Psi,Psi2);
// Psi = Give_coord_quadrangle(quadIter,quadAddX,quadAddY+i*quadFix[1]);
}
if (id == quad->Give_id_edge(Edir2D))
{
// cout << "Quadrangle id: " << quadIter << " ";
// cout << "correct edge found: Edir2D " << endl;
// Give_coord_quadrangle(quadIter,i,0).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,0,i).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 )-i,0).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,0,Give_Nx_hexahedron ( 0 )-i).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,i,Give_Nx_hexahedron ( 0 )).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 ),i).Print();cout<<endl; // korrekt
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 )-i,Give_Nx_hexahedron ( 0 )).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 ),Give_Nx_hexahedron ( 0 )-i).Print();cout<<endl;
// cout << endl;
Psi = Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 ),i);
//D3vector Psi2 = Give_coord_quadrangle(quadIter,Nx-quadAddX,quadAddY+i*quadFix[1]);
//compareVectors(Psi,Psi2);
//Psi = Give_coord_quadrangle(quadIter,Nx-quadAddX,quadAddY+i*quadFix[1]);
}
if (id == quad->Give_id_edge(Sdir2D))
{
// cout << "Quadrangle id: " << quadIter << " ";
// cout << "correct edge found: Sdir2D " << endl;
// Give_coord_quadrangle(quadIter,i,0).Print();cout<<endl; //korrekt
// Give_coord_quadrangle(quadIter,0,i).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 )-i,0).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,0,Give_Nx_hexahedron ( 0 )-i).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,i,Give_Nx_hexahedron ( 0 )).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 ),i).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 )-i,Give_Nx_hexahedron ( 0 )).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 ),Give_Nx_hexahedron ( 0 )-i).Print();cout<<endl;
Psi = Give_coord_quadrangle(quadIter,i,0);
//D3vector Psi2 = Give_coord_quadrangle(quadIter,quadAddX+i*quadFix[0],quadAddY);
//compareVectors(Psi,Psi2);
//Psi = Give_coord_quadrangle(quadIter,quadAddX+i*quadFix[0],quadAddY);
}
if (id == quad->Give_id_edge(Ndir2D))
{
// cout << "Quadrangle id: " << quadIter << " ";
// cout << "correct edge found: Ndir2D " << endl;
// Give_coord_quadrangle(quadIter,i,0).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,0,i).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 )-i,0).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,0,Give_Nx_hexahedron ( 0 )-i).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,i,Give_Nx_hexahedron ( 0 )).Print();cout<<endl; // KORREKT
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 ),i).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 )-i,Give_Nx_hexahedron ( 0 )).Print();cout<<endl;
// Give_coord_quadrangle(quadIter,Give_Nx_hexahedron ( 0 ),Give_Nx_hexahedron ( 0 )-i).Print();cout<<endl;
Psi = Give_coord_quadrangle(quadIter,i,Give_Nx_hexahedron ( 0 ));
// D3vector Psi2 = Give_coord_quadrangle(quadIter,quadAddX+i*quadFix[0],Ny - quadAddY);
//compareVectors(Psi,Psi2);
//Psi = Give_coord_quadrangle(quadIter,quadAddX+i*quadFix[0],Ny - quadAddY);
}
}
return Psi;
cout << endl;
return Psi;
// cout << endl;
//return
......@@ -1067,7 +1177,7 @@ D3vector Blockgrid::Give_coord_edge ( int id, int i ) const
cout << "number quadra : "<<ug->Give_number_quadrangles()<<endl;
//Christoph version
}
double eta;
D3vector PL, PR;
Edge_el* edge;
......@@ -1076,7 +1186,7 @@ D3vector Blockgrid::Give_coord_edge ( int id, int i ) const
eta = ( double ) i/ ( double ) Give_Nx_edge ( id );
cout << "psi new "; Psi.Print();cout << endl;
//cout << "psi new "; Psi.Print();cout << endl;
Psi = edge->transform ( eta , edge->shiftPointer(ug->Get_pointer_global_data()));
......@@ -1088,7 +1198,7 @@ D3vector Blockgrid::Give_coord_edge ( int id, int i ) const
PR = ug->Give_point(edge->Give_id_corner_E())->Give_coordinate();
*/
cout << "psi old " ; (Psi + PL + ( PR-PL ) * eta).Print();cout << endl;
// cout << "psi old " ; (Psi + PL + ( PR-PL ) * eta).Print();cout << endl;
return Psi + PL + ( PR-PL ) * eta;
}
......
......@@ -104,6 +104,7 @@ class Blockgrid {
private:
Variable<std::complex<double> > *phase_shift;
bool transformFromQuadrangle;
bool variable_set; // set true falls eine Variable konstruiert
int id_of_grid;
......
......@@ -400,11 +400,7 @@ D3vector transform_left_lens_diag_NE_quad ( double t1, double t2, double* pointe
double radiusSquared = ((x+actualX)*(x+actualX)+(y+actualY)*(y+actualY));
double z = offsetZ_global_data+
sign*( ( sqrt(pow(curvatureRight_global_data,2)-radiusSquared)
+ sign*(-curvatureRight_global_data +( (1-t)*(thickness_global_data-z_right_inner_global_data) + t * (thickness_global_data-z_right_outer_global_data)))));
z =offsetZ_global_data+sign*(-1*( sqrt(pow(curvatureLeft_global_data,2)-radiusSquared) +sign *(-curvatureLeft_global_data +( (1-t)*z_left_inner_global_data + t * z_left_outer_global_data))));
double z =offsetZ_global_data+sign*(-1*( sqrt(pow(curvatureLeft_global_data,2)-radiusSquared) +sign *(-curvatureLeft_global_data +( (1-t)*z_left_inner_global_data + t * z_left_outer_global_data))));
if (std::isnan(z))
{z = 0;}
......@@ -422,7 +418,7 @@ D3vector transform_diag_inner_faces_NE_quad( double t1, double t2, double* point
double z_left_outer_global_data = pointer_global_data[7];
double z_right_inner_global_data = pointer_global_data[8];
double z_right_outer_global_data = pointer_global_data[9];
bool invertT1T2 = false;
bool invertT1T2 = true;
if (invertT1T2)
{
......@@ -446,9 +442,14 @@ D3vector transform_diag_inner_faces_NE_quad( double t1, double t2, double* point
sign*( ( sqrt(pow(curvatureRight_global_data,2)-
( radiusSquared))
+ sign*(-curvatureRight_global_data +( (1-t)*(thickness_global_data-z_right_inner_global_data) + t * (thickness_global_data-z_right_outer_global_data)))));
double zLeft = offsetZ_global_data
+sign*(-1*( sqrt(pow(curvatureLeft_global_data,2)-
radiusSquared)
+sign *(-curvatureLeft_global_data +( (1-t)*z_left_inner_global_data + t * z_left_outer_global_data))));
double z = zRight * (1-t1) ;
z = 0;
double z = zRight * (t1) + zLeft * (1-t1) ;
if (( t1 == 0 && t2 == 0)||( t1 == 1 && t2 == 0)||( t1 == 0 && t2 == 1)||( t1 == 1 && t2 == 1))
{
if ( fabs(z) > 1e-10)
......
......@@ -947,6 +947,27 @@ void Unstructured_grid::construction_done() {
// IV. Zordering
zordering = new Zordering(this);
//zordering->Print(this);
// V. Check, whether FACE or EDGE transform is used to find coordinates. Default is FALSE
transformFromQuadrangle = true;
//transformFromQuadrangle
for (int qu = 0; qu < Give_number_quadrangles(); ++qu)
{
if( Give_quadrangle(qu)->transform == NULL)
{
transformFromQuadrangle = false;
}
}
if (transformFromQuadrangle)
{ edge_to_quad_id.resize(Give_number_edges());
edge_to_quad_dir.resize(Give_number_edges());
for (int iter = 0 ; iter <edge_to_quad_id.size(); iter++ )
{
edge_to_quad_id.at(iter) = -1;
}}
}
void Unstructured_grid::Set_recursive_color(int id_edge, int color) {
......
......@@ -110,7 +110,6 @@ class Zordering {
private:
int* indices;
elementTyp* type_index;
int number_indizes;
};
......@@ -168,8 +167,15 @@ class Unstructured_grid : public Partitioning {
int Give_number_edges() { return num_edges; }
int Give_number_quadrangles() { return num_quadrangles; }
bool Give_transform_From_Quadrangle() {return transformFromQuadrangle;}
// bool Adjacent(int num_hex, int num_edge);
void Set_edge_to_quad_id(int index, int value){ if(index < Give_number_edges());edge_to_quad_id.at(index) = value; };
int Give_edge_to_quad_id(int index){ if(index < Give_number_edges());return edge_to_quad_id.at(index); };
void Set_edge_to_quad_dir(int index, int value){ if(index < Give_number_edges());edge_to_quad_dir.at(index) = value; };
int Give_edge_to_quad_dir(int index){ if(index < Give_number_edges());return edge_to_quad_dir.at(index); };
Hexahedron_el* Give_hexahedron(int i) { return &(hexahedra[i]); }
Quadrangle_el* Give_quadrangle(int i) { return &(quadrangles[i]); }
Edge_el* Give_edge(int i) { return &(edges[i]); }
......@@ -247,6 +253,7 @@ class Unstructured_grid : public Partitioning {
std::vector<double> construction_lx, construction_ly, construction_lz;
bool constructionBoolArched;
bool constructionBoolPeriodic;
bool transformFromQuadrangle;
private:
......@@ -254,7 +261,8 @@ class Unstructured_grid : public Partitioning {
bool construction_done_yn;
bool periodic; // is num_poi == num_poi_plus_periodic
std::vector<int> edge_to_quad_id;
std::vector<int> edge_to_quad_dir;
Topf<2, sorted_pair> topf_edges; // spaeter löschen?
Topf<4, sorted_quad> topf_quadrangles; // spaeter löschen?
......
Markdown is supported
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