continuous Arg calculator eingebaut

program2D/source/extemp/calcArg.cc

+/**********************************************************************************
+ * Copyright 2010 Christoph Pflaum 
+ * 		Department Informatik Lehrstuhl 10 - Systemsimulation
+ *		Friedrich-Alexander Universität Erlangen-Nürnberg
+ * 
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ **********************************************************************************/
+ 
+ 
+// ------------------------------------------------------------
+//
+// calcArg.h
+//
+// ------------------------------------------------------------
+
+#include <string>
+#include <stdexcept>
+#include <math.h>
+
+#include "calcArg.h"
+#include "../abbrevi2D.h"
+
+#include "../math_lib/math_lib.h"
+#include "../grid/elements2D.h"
+#include "../grid/parti2D.h"
+#include "../grid/ug2D.h"
+#include "../grid/blockgrid2D.h"
+#include "../grid/marker2D.h"
+#include "extemp2D.h"
+#include "variable2D.h"
+#include "variable2D_cc.h"
+#include "update2D.h"
+#include "print_var2D_vtk_cc.h"
+#include "co_fu2D.h"
+
+
+PointForPhase::PointForPhase(int i_, int j_, double r_, int iSmall_, int jSmall_) {
+    i = i_;
+    j = j_;
+    r = r_;
+    iSmall = iSmall_;
+    jSmall = jSmall_;
+}
+
+bool compare_Radius(PointForPhase& pA, PointForPhase& pB) {
+   return  pA.getRadius() < pB.getRadius();
+}
+
+bool compare_Radius_blocks(BlockForPhase& pA, BlockForPhase& pB) {
+   return  pA.getRadius() < pB.getRadius();
+}
+
+CalcContinuousArg::CalcContinuousArg(Blockgrid2D& blockgrid_) {  
+    blockgrid = &blockgrid_;
+    ug = blockgrid->Give_unstructured_grid();
+    
+    Variable2D<double> norm(blockgrid_);
+    X_coordinate2d X(blockgrid_);
+    Y_coordinate2d Y(blockgrid_);
+    
+    norm = X*X+Y*Y;
+    
+    numberBlocks = ug->Give_number_rectangle();
+    
+    orderedPoints = new std::list<PointForPhase>[numberBlocks];
+        
+    enum woSmaller { nichts, Norden, Sueden, Osten, Westen };
+        
+    for(int id=0;id<numberBlocks;++id) {
+        int Nx = blockgrid->Give_Nx_rectangle(id);
+        int Ny = blockgrid->Give_Ny_rectangle(id);
+        for(int i=0;i<=Nx;++i) for(int j=0;j<=Ny;++j) { 
+	    norm.Update<rectangleEl>(id);
+	    woSmaller woNext = nichts;
+	    double radiusNext = -1.0;
+	    double radiusM = norm.template Give_data<rectangleEl>(id,i,j,Nx);
+	    
+	    //Westen
+            if(i>0) { 
+	       double r = norm.template Give_data<rectangleEl>(id,i-1,j,Nx);
+               if(r < radiusM && r > radiusNext) {
+		  radiusNext = r;
+		  woNext = Westen;
+	       }
+	    }
+	    //Osten
+            if(i<Nx) { 
+	       double r = norm.template Give_data<rectangleEl>(id,i+1,j,Nx);
+               if(r < radiusM && r > radiusNext) {
+		  radiusNext = r;
+		  woNext = Osten;
+	       }
+	    }
+	    //Sueden
+            if(j>0) { 
+	       double r = norm.template Give_data<rectangleEl>(id,i,j-1,Nx);
+               if(r < radiusM && r > radiusNext) {
+		  radiusNext = r;
+		  woNext = Sueden;
+	       }
+	    }
+	    //Norden
+            if(i<Ny) { 
+	       double r = norm.template Give_data<rectangleEl>(id,i,j+1,Nx);
+               if(r < radiusM && r > radiusNext) {
+		  radiusNext = r;
+		  woNext = Norden;
+	       }
+	    }
+	    if(woNext == Westen) orderedPoints[id].push_back(PointForPhase(i,j,radiusM,i-1,j));
+	    if(woNext == Osten)  orderedPoints[id].push_back(PointForPhase(i,j,radiusM,i+1,j));
+	    if(woNext == Sueden) orderedPoints[id].push_back(PointForPhase(i,j,radiusM,i,j-1));
+	    if(woNext == Norden) orderedPoints[id].push_back(PointForPhase(i,j,radiusM,i,j+1));   
+	}
+	orderedPoints[id].sort(compare_Radius);    
+	
+	orderedBlocks.push_back(BlockForPhase(id,orderedPoints[id].front().getRadius()));
+    }
+    
+    orderedBlocks.sort(compare_Radius_blocks);    
+    
+     for(std::list<BlockForPhase>::iterator it=orderedBlocks.begin(); it != orderedBlocks.end(); ++it) {
+         std::cout << " id: " << it->getId() << " R: " << it->getRadius() << std::endl;
+     }
+}
+
+double phaseComp(std::complex<double> value) {
+  return atan2(value.imag(),value.real());
+}
+
+void CalcContinuousArg::calcArg(Variable2D<double>& phase, Variable2D<std::complex<double> >& Efield) {
+   Function2d1<double, std::complex<double> >  Phase(phaseComp);
+   phase = Phase(Efield);           
+  
+   for(std::list<BlockForPhase>::iterator it=orderedBlocks.begin(); it != orderedBlocks.end(); ++it) {
+       int id = it->getId();
+       
+       for(std::list<PointForPhase>::iterator itP=orderedPoints[id].begin(); itP != orderedPoints[id].end(); ++itP) {
+	   phase.Update<rectangleEl>(id);
+           int Nx = blockgrid->Give_Nx_rectangle(id);       
+	 
+	   int i = itP->getI();
+	   int j = itP->getJ();
+	   
+	   std::complex<double> EField_M = Efield.template Give_data<rectangleEl>(id,i,j,Nx);
+            
+           double phaseM = atan2(EField_M.imag(),EField_M.real());
+	   
+	   int iSmall = itP->getIsmall();
+	   int jSmall = itP->getJsmall();
+           
+           double phaseNext = phase.template Give_data<rectangleEl>(id,iSmall,jSmall,Nx);
+ 	   
+	   if(fabs(phaseM-phaseNext)>2.0) {
+	      if(phaseM < phaseNext) while(phaseM < phaseNext) phaseM = phaseM + 2.0*M_PI;
+	      else                   while(phaseM > phaseNext) phaseM = phaseM - 2.0*M_PI;
+              phase.setValueAtRectangularPoint(phaseM,id,i,j,Nx);
+	   }	   
+       }
+       phase.Update_back(id);
+    }
+}
+
+
+
+
+
+

program2D/source/extemp/calcArg.h

+/**********************************************************************************
+ * Copyright 2010 Christoph Pflaum 
+ * 		Department Informatik Lehrstuhl 10 - Systemsimulation
+ *		Friedrich-Alexander Universität Erlangen-Nürnberg
+ * 
+ * Licensed under the Apache License, Version 2.0 (the "License");
+ * you may not use this file except in compliance with the License.
+ * You may obtain a copy of the License at
+ *
+ * http://www.apache.org/licenses/LICENSE-2.0
+ *
+ * Unless required by applicable law or agreed to in writing, software
+ * distributed under the License is distributed on an "AS IS" BASIS,
+ * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+ * See the License for the specific language governing permissions and
+ * limitations under the License.
+ **********************************************************************************/
+ 
+// ------------------------------------------------------------
+//
+// calcArg.h
+//
+// ------------------------------------------------------------
+#include "../grid/ug2D.h"
+#include "../grid/blockgrid2D.h"
+
+#ifndef CALCARG2D_H
+#define CALCARG2D_H
+#ifdef _OPENMP
+#include <omp.h>
+#include "../ugopenmp.h"
+#endif
+
+#include <list>
+
+class PointForPhase {
+  public:
+    PointForPhase(int i_, int j_, double r_, int iSmall_, int jSmall_);
+    double getRadius() { return r; }
+    
+    int getI() { return i; }
+    int getJ() { return j; }
+    int getIsmall() { return iSmall; }
+    int getJsmall() { return jSmall; }    
+  private:
+    int i;
+    int j;
+    double r;
+    int iSmall;
+    int jSmall;
+};
+
+
+class BlockForPhase {
+  public:
+    BlockForPhase(int id_, double r_) { id = id_;  r = r_; }
+    double getRadius() { return r; }
+    double getId() { return id; }      
+  private:
+    int id;
+    double r;
+};
+
+class CalcContinuousArg {
+  public:  
+    CalcContinuousArg(Blockgrid2D& blockgrid_);
+    void calcArg(Variable2D<double>& phase, Variable2D<std::complex<double> >& Efield);    
+    
+  private:
+    Blockgrid2D*        blockgrid;
+    Unstructured2DGrid* ug;
+    std::list<PointForPhase>* orderedPoints; 
+    std::list<BlockForPhase> orderedBlocks; 
+    int numberBlocks;
+   
+    
+};
+
+
+
+#endif // CALCARG2D_H
+
+
+

program2D/source/extemp/variable2D.h

@@ -159,10 +159,13 @@ template <class DTyp>                       class Variable2D;
 template <class A, class DTyp, class DTyp1> class Exp_Function1;
 template <class A, class B, class C, class D, class DTyp, class DTyp_Result> class Exp_Function4;
 
+class CalcContinuousArg;
+
 /** \addtogroup ExpressionTemplates2D **/
 /* @{ */  
 template <class DTyp>
 class Variable2D : public Expr2D<Variable2D< DTyp > > {
+        friend class CalcContinuousArg;
         public:
                 typedef DTyp Result;
 
@@ -202,7 +205,7 @@ class Variable2D : public Expr2D<Variable2D< DTyp > > {
 
                 template <elementTyp TYP_EL>
                 inline DTyp Give_data ( params2D_in ) const;
-
+		
 		// fuer debug-zwecke mit positionsrueckgabe
                 template <elementTyp TYP_EL>
                 inline DTyp Give_data ( params2D_in, D2vector &posVector ) const ;
@@ -268,6 +271,8 @@ class Variable2D : public Expr2D<Variable2D< DTyp > > {
                 //template <elementTyp TYP_EL>
                 void Update_back(int id) const;
 
+		void setValueAtRectangularPoint(DTyp val, int id, int i, int j, int Nx) { data_rectangles[id][i+ ( Nx+1 ) *j] = val; }		
+		
                 Blockgrid2D*        blockgrid;
                 Unstructured2DGrid*        ug;
 

program2D/source/ugblock2D.h

@@ -43,6 +43,7 @@
 
 
 #include "extemp/extemp2D.h"
+#include "extemp/calcArg.h"
 #include "extemp/co_fu2D.h"
 #include "extemp/functor2D.h"
 #include "extemp/variable2D.h"

program2D/ugblocks2D.pri

@@ -28,6 +28,7 @@ INCLUDEPATH += . \
 # Input
 HEADERS += $${UGBLOCKS_PATH}/program2D/source/abbrevi2D.h \
            $${UGBLOCKS_PATH}/program2D/source/ugblock2D.h \
+           $${UGBLOCKS_PATH}/program2D/source/extemp/calcArg.h \           
            $${UGBLOCKS_PATH}/program2D/source/extemp/cellvar2D.h \
            $${UGBLOCKS_PATH}/program2D/source/extemp/co_fu2D.h \
            $${UGBLOCKS_PATH}/program2D/source/extemp/functor2D.h \
@@ -61,6 +62,7 @@ HEADERS += $${UGBLOCKS_PATH}/program2D/source/abbrevi2D.h \
   
 SOURCES += $${UGBLOCKS_PATH}/program2D/source/extemp/variable2D.cc \
            $${UGBLOCKS_PATH}/program2D/source/extemp/cellvar2D.cc \
+           $${UGBLOCKS_PATH}/program2D/source/extemp/calcArg.cc \                      
            $${UGBLOCKS_PATH}/program2D/source/extemp/print_var2D_Qt_vtk.cc \
            $${UGBLOCKS_PATH}/program2D/source/extemp/extemp2D.cc \           
            $${UGBLOCKS_PATH}/program2D/source/grid/blockgrid2D.cc \