src/expde/grid/p_opendx.cc

Go to the documentation of this file.

//    expde: expression templates for partial differential equations.
//    Copyright (C) 2001  Christoph Pflaum
//    This program is free software; you can redistribute it and/or modify
//    it under the terms of the GNU General Public License as published by
//    the Free Software Foundation; either version 2 of the License, or
//    (at your option) any later version.
//
//    This program is distributed in the hope that it will be useful,
//    but WITHOUT ANY WARRANTY; without even the implied warranty of
//    MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
//    GNU General Public License for more details.
//
//                 SEE  Notice1.doc made by 
//                 LAWRENCE LIVERMORE NATIONAL LABORATORY
//

// ------------------------------------------------------------
//
// p_OpenDx.cc
//
// ------------------------------------------------------------

// for OpenDx


// Include level 0:
#ifdef COMP_GNUOLD
 #include <iostream.h>
 #include <fstream.h>
 #include <time.h>
 #include <math.h>
#else
 #include <iostream>
 #include <fstream>
 #include <ctime>
 #include <cmath>
#endif 

  
// Include level 0:
#include "../parser.h"

// Include level 1:
#include "../paramete.h"
#include "../abbrevi.h"
#include "../math_lib/math_lib.h"
 
// Include level 2:
#include "../basic/basic.h"

// Include level 3:
#include "../domain/domain.h"

// Include level 4:
#include "../formulas/boundy.h"
#include "../formulas/loc_sten.h"

// Include level 5:
#include "gpar.h"
#include "parallel.h"
#include "mgcoeff.h"
#include "sto_man.h"
#include "gridbase.h"
#include "grid.h"

// a)  for AVS

void Grid_base::Print_ec_opendx(IndexBo indexbo_A, IndexBo indexbo_B,
                                   ofstream *Datei, int number) {
  *Datei << number << "  1  line  ";
  *Datei << hashtable3_point(indexbo_A)->nummer  << "  ";
  *Datei << hashtable3_point(indexbo_B)->nummer  << "  ";
  *Datei << endl; 
}

// Ausgabe der Tetraeder
int Grid_base::opendx_bo_cell(BoCell* bo,ofstream *datei,bool print_or_calc, 
                           int number) {
  Tetraeder_storage* tets;
  int num;
  Index3D I;
  Cell_type_points typ;

#ifdef PERIODIC
  if(bo->Give_Index().neighbour_non_periodic(WSDd).Is_non_periodic()) {
    return number;
  }
#endif


  if(print_or_calc==false) {
    for(tets = bo->Give_tets();tets!=NULL;tets = tets->Next()) {
      number++;
    }
  }
  else {
    I = bo->Give_Index();

    for(tets = bo->Give_tets();tets!=NULL;tets = tets->Next()) {
      ++number; 
      // Punkt 0:
      num = tets->N0();
      typ = bo->edge_point(num);
      if(typ==edge_poi_typ) {
        *datei << Give_Nummer(I.neighbour(bo->corner(num)),
                              bo->edge_dir(num)) << "  ";
      }
      if(typ==corner_poi_typ) {
        *datei << Give_Nummer(I.neighbour(bo->corner(num)))
               << "  "; 
      }
      if(typ==cell_poi_typ) {
        *datei << Give_Nummer_cellpoi(I) << "  ";       
      }
      // Punkt 1:
      num = tets->N1();
      typ = bo->edge_point(num);
      if(typ==edge_poi_typ) {
        *datei << Give_Nummer(I.neighbour(bo->corner(num)),
                              bo->edge_dir(num)) << "  ";
      }
      if(typ==corner_poi_typ) {
        *datei << Give_Nummer(I.neighbour(bo->corner(num)))
               << "  "; 
      }
      if(typ==cell_poi_typ) {
        *datei << Give_Nummer_cellpoi(I) << "  ";       
      }
      // Punkt 2:
      num = tets->N2();
      typ = bo->edge_point(num);
      if(typ==edge_poi_typ) {
        *datei << Give_Nummer(I.neighbour(bo->corner(num)),
                              bo->edge_dir(num)) << "  ";
      }
      if(typ==corner_poi_typ) {
        *datei << Give_Nummer(I.neighbour(bo->corner(num)))
               << "  "; 
      }
      if(typ==cell_poi_typ) {
        *datei << Give_Nummer_cellpoi(I) << "  ";       
      }
      // Punkt 3:
      num = tets->N3();
      typ = bo->edge_point(num);
      if(typ==edge_poi_typ) {
        *datei << Give_Nummer(I.neighbour(bo->corner(num)),
                              bo->edge_dir(num)) << "  ";
      }
      if(typ==corner_poi_typ) {
        *datei << Give_Nummer(I.neighbour(bo->corner(num)))
               << "  "; 
      }
      if(typ==cell_poi_typ) {
        *datei << Give_Nummer_cellpoi(I) << "  ";       
      }
      *datei << endl;
    }
  }
  return number;
}



int Grid_base::Recursion_Cells_OPENDX(Index3D I, ofstream *Datei, int nummer) {
  int i;
  Index3D Ison;
  for(i=0;i<8;++i) {
    if(Exists_Point(I.son((dir_sons)i)))
      nummer = Recursion_Cells_OPENDX(I.son((dir_sons)i),Datei,nummer);
    else {
      Ison = I.son((dir_sons)i);

     // periodic
#ifdef PERIODIC
      if(Ison.neighbour_non_periodic(WSDd).Is_non_periodic()
         && Give_cell_typ(Ison) == int_cell) {
#else
      if(Give_cell_typ(Ison) == int_cell) {
#endif
        nummer++;

        // 0. Tet: WND, WNT, WST, EST 
        *Datei << Give_Nummer(Ison.neighbour(WNTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(WNDd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(WSTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(ESTd))  << "  ";
        *Datei << "\n";
        
        nummer++;
        // 1. Tet: EST, WND, WST, ESD  
        *Datei << Give_Nummer(Ison.neighbour(ESTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(WNDd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(WSTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(ESDd))  << "  ";
        *Datei << "\n";
        
        nummer++;
        // 2. Tet: WND, WSD, WST, ESD        
        *Datei << Give_Nummer(Ison.neighbour(WNDd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(WSDd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(WSTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(ESDd))  << "  ";
        *Datei << "\n";
        
        nummer++;
        // 3. Tet: EST, WND, ESD, END
        *Datei << Give_Nummer(Ison.neighbour(ESTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(WNDd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(ESDd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(ENDd))  << "  ";
        *Datei << "\n";
        
        nummer++;
        // 4. Tet: ENT, WNT, EST, END
        *Datei << Give_Nummer(Ison.neighbour(ENTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(WNTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(ESTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(ENDd))  << "  ";
        *Datei << "\n";
        
        nummer++;
        // 5. Tet: WNT, WND, EST, END
        *Datei << Give_Nummer(Ison.neighbour(WNTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(WNDd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(ESTd))  << "  ";
        *Datei << Give_Nummer(Ison.neighbour(ENDd))  << "  ";
        *Datei << "\n";
      }
    }
  }
  return nummer;
}


// There are the main print functions
void Grid_base::Print_Variable_OpenDx(ofstream *Datei, int number_var) {
  int number_cells, number, ebene;
  int number_points;
  Index3D ind;
  Pointtype typ_point;

  if(Is_storage_initialized()==false) { 
    cout << " \n Fehler in Print_Variable_OpenDx! Initialisierung fehlt! "
         << endl;
    return;
  }
  if(Give_max_num_var() <= number_var) {
    cout << " \n Fehler in Print_Variable_OPENDX! number_var zu gross! " << endl;
  }
  else {
    // Teil 0: Write information
    *Datei << "# File created by ExPDE " << endl;
    *Datei << "# dx file format for OpenDx " << endl;
    *Datei << "# scalar value " << endl;

    // Teil 1: Anzahl der Zellen zaehlen
    // Teil 1.1: Innerer Zellen
    ind = Index3D(2,2,2);
    number_cells = Recursion_Count_Cells(ind);
    // Teil 1.2: Randzellen
    iterate_hash2 {
      number_cells = opendx_bo_cell(bocell,Datei,false,number_cells);
    }
    // Teil 1.3: Anzahl der Punkte
    // a) innen und randnah
    number_points=0;
    iterate_hash1 {
      if(point1->typ >= interior) {
        point1->nummer = number_points;
        number_points++;
      }
      else point1->nummer = -1;
    }
    // b) Randpunkte
    iterate_hash3 {
      bo2point->Set_number_avs(number_points);
      ++number_points;
    }
    // c) Randzellfreiheitsgrad
    iterate_hash2 {
      if(bocell->Exists_bocellpoint()) {
        bocell->Set_number_avs(number_points);
        number_points++;
      }
    }

    // Teil 2: Koordinaten der Punkte eingeben
    // Teil 2.0: Header file for points
    *Datei << "object 1 class array type float rank 1 shape 3 items "
           << number_points << " data follows" << endl;

    // Teil 2.1 alle inneren Punkte
    iterate_hash1 {
      if(point1->nummer != -1) {
        *Datei << " " ;
        transform_coord(point1->ind.coordinate()).Print(Datei);
        *Datei  << "\n";
      }
    }
    // Teil 2.2 alle Randpunkte (Kantenpunkte)
    iterate_hash3 {
      *Datei << " " ;
      (bo2point->transform_coord(Give_A(),
                                 H_mesh())).Print(Datei);
      *Datei  << "\n";
    }
    // Teil 2.3 alle Randzellfreiheitsgrad
    iterate_hash2 {
      if(bocell->Exists_bocellpoint()) {
        *Datei << " " ;
        (transform_coord(bocell->ind.coordinate()) + 
         bocell->Local_coord_bocellpoint()).Print(Datei);
        *Datei  << "\n";
      }
    }
    // Teil 3: Zellen ausgeben
    // Teil 3.0: Header file for tetrahedra:
    *Datei << "object 2 class array type int rank 1 shape 4 items "
           << number_cells << " data follows" << endl;
    // Teil 3.1: Innere Zellen ausgeben
    ind = Index3D(2,2,2);
    number=Recursion_Cells_OPENDX(ind,Datei,0);
    // Teil 3.2: Randzellen
    iterate_hash2 {
      number = opendx_bo_cell(bocell,Datei,true,number);
    }
    // Teil 3.3: end for tetrahedra
    *Datei << "attribute " << '"' << "element type" 
           << '"' << " string " << '"' << "tetrahedra" << '"' << "" << endl;
    *Datei << "attribute " << '"' << "ref" << '"' 
           << " string " << '"' << "positions" << '"' << "" << endl;

    // Teil 4: Header file for one variable:
    *Datei << "object 3 class array type float rank 0 items "
           << number_points << " data follows" << endl;
    
    // Teil 6: Ausgabe der Werte der Variablen
    // Teil 6.1: Werte im Inneren
    iterate_hash1 {
      if(point1->nummer != -1) {
        //      *Datei << point1->nummer << "  ";
        typ_point = point1->typ;
        if(typ_point==exterior)
          cout << " Fehler in Print_Variable_OPENDX! " << endl;
        if(typ_point!=multigrid) {
          ebene = Give_finest_level(point1->ind);
          *Datei << Give_variable(point1->ind,ebene)[number_var];
          *Datei << "\n";
        }
        else *Datei << " 1.0 \n";
      }
    }
    // Teil 6.2: Werte am Rand
    iterate_hash3 {
      //  *Datei << bo2point->nummer+hashtable1_occ << "  ";
      *Datei << Give_variable(bo2point->ind,bo2point->direction)[number_var];
      *Datei << " \n";
    }
    // Teil 6.3 Werte bei Randzellfreiheitsgrad
    iterate_hash2 {
      if(bocell->Exists_bocellpoint()) {
        //      *Datei << bocell->Number_avs() << " " ;
        *Datei << bocell->Give_var()[number_var];
        *Datei  << "\n";
      }
    }
    // Teil 6.4: end variable
    *Datei << "attribute " 
           << '"' << "dep" 
           << '"' << " string " << '"' << "positions" << '"' << "" << endl;

    // Teil 7: end remarks
    *Datei << "object " << '"' 
           << "irregular positions irregular connections"
           << '"' << " class field" << 
      endl;
    *Datei << "component " << '"' << "positions" << '"' << " value 1" << endl;
    *Datei << "component " 
           << '"' << "connections" << '"' << " value 2" << endl;
    *Datei << "component " << '"' << "data" << '"' << " value 3" << endl;
    *Datei << "end" << endl; 

//cout << " Hi: no point "; 

  }
}




void Grid_base::Print_Variable_OpenDx_moved(ofstream *Datei, int number_var,
                                            int number_var_a,
                                            int number_var_b, 
                                            int number_var_c) {
  int number_cells, number, ebene;
  int number_points;
  Index3D ind;
  Pointtype typ_point;

  if(Is_storage_initialized()==false) { 
    cout << 
      " \n Fehler in Print_Variable_OpenDx_moved! Initialisierung fehlt! "
         << endl;
    return;
  }
  if(Give_max_num_var() <= number_var) {
    cout << " \n Fehler in Print_Variable_OPENDX! number_var zu gross! " << endl;
  }
  else {
    // Teil 0: Write information
    *Datei << "# File created by ExPDE " << endl;
    *Datei << "# dx file format for OpenDx " << endl;
    *Datei << "# scalar value moved by a vector" << endl;

    // Teil 1: Anzahl der Zellen zaehlen
    // Teil 1.1: Innerer Zellen
    ind = Index3D(2,2,2);
    number_cells = Recursion_Count_Cells(ind);
    // Teil 1.2: Randzellen
    iterate_hash2 {
      number_cells = opendx_bo_cell(bocell,Datei,false,number_cells);
    }
    // Teil 1.3: Anzahl der Punkte
    // a) innen und randnah
    number_points=0;
    iterate_hash1 {
      if(point1->typ >= interior) {
        point1->nummer = number_points;
        number_points++;
      }
      else point1->nummer = -1;
    }
    // b) Randpunkte
    iterate_hash3 {
      bo2point->Set_number_avs(number_points);
      ++number_points;
    }
    // c) Randzellfreiheitsgrad
    iterate_hash2 {
      if(bocell->Exists_bocellpoint()) {
        bocell->Set_number_avs(number_points);
        number_points++;
      }
    }
    // Teil 2: Koordinaten der Punkte eingeben
    // Teil 2.0: Header file for points
    *Datei << "object 1 class array type float rank 1 shape 3 items "
           << number_points << " data follows" << endl;

    // Teil 2.1 alle inneren Punkte
    iterate_hash1 {
      if(point1->nummer != -1) {
        *Datei << " " ;
        (transform_coord(point1->ind.coordinate()) +
         D3vector(Give_variable(point1->ind,Max_level())[number_var_a],
                  Give_variable(point1->ind,Max_level())[number_var_b],
                  Give_variable(point1->ind,Max_level())[number_var_c])
         ).Print(Datei);
        *Datei  << "\n";
      }
    }
    // Teil 2.2 alle Randpunkte (Kantenpunkte)
    iterate_hash3 {
      *Datei << " " ;
      ((bo2point->transform_coord(Give_A(),
                                  H_mesh())) +
       D3vector(Give_variable(bo2point->ind,bo2point->direction)[number_var_a],
                Give_variable(bo2point->ind,bo2point->direction)[number_var_b],
                Give_variable(bo2point->ind,bo2point->direction)[number_var_c])
       ).Print(Datei);
      *Datei  << "\n";
    }
    // Teil 2.3 alle Randzellfreiheitsgrad
    iterate_hash2 {
      if(bocell->Exists_bocellpoint()) {
        *Datei << " " ;
        (transform_coord(bocell->ind.coordinate()) + 
         bocell->Local_coord_bocellpoint() +
         D3vector(bocell->Give_var()[number_var_a],
                  bocell->Give_var()[number_var_b],
                  bocell->Give_var()[number_var_c])
         ).Print(Datei);
        *Datei  << "\n";
      }
    }
    // Teil 3: Zellen ausgeben
    // Teil 3.0: Header file for tetrahedra:
    *Datei << "object 2 class array type int rank 1 shape 4 items "
           << number_cells << " data follows" << endl;
    // Teil 3.1: Innere Zellen ausgeben
    ind = Index3D(2,2,2);
    number=Recursion_Cells_OPENDX(ind,Datei,0);
    // Teil 3.2: Randzellen
    iterate_hash2 {
      number = opendx_bo_cell(bocell,Datei,true,number);
    }
    // Teil 3.3: end for tetrahedra
    *Datei << "attribute " << '"' << "element type" 
           << '"' << " string " << '"' << "tetrahedra" << '"' << "" << endl;
    *Datei << "attribute " << '"' << "ref" << '"' 
           << " string " << '"' << "positions" << '"' << "" << endl;

    // Teil 4: Header file for one variable:
    *Datei << "object 3 class array type float rank 0 items "
           << number_points << " data follows" << endl;
    
    // Teil 6: Ausgabe der Werte der Variablen
    // Teil 6.1: Werte im Inneren
    iterate_hash1 {
      if(point1->nummer != -1) {
        //      *Datei << point1->nummer << "  ";
        typ_point = point1->typ;
        if(typ_point==exterior)
          cout << " Fehler in Print_Variable_OPENDX! " << endl;
        if(typ_point!=multigrid) {
          ebene = Give_finest_level(point1->ind);
          *Datei << Give_variable(point1->ind,ebene)[number_var];
          //     *Datei << 2.0;
          *Datei << "\n";
        }
        else *Datei << " 1.0 \n";
      }
    }
    // Teil 6.2: Werte am Rand
    iterate_hash3 {
      //  *Datei << bo2point->nummer+hashtable1_occ << "  ";
      *Datei << Give_variable(bo2point->ind,bo2point->direction)[number_var];
      *Datei << " \n";
    }
    // Teil 6.3 Werte bei Randzellfreiheitsgrad
    iterate_hash2 {
      if(bocell->Exists_bocellpoint()) {
        //      *Datei << bocell->Number_avs() << " " ;
        *Datei << bocell->Give_var()[number_var];
        *Datei  << "\n";
      }
    }
    // Teil 6.4: end variable
    *Datei << "attribute " 
           << '"' << "dep" 
           << '"' << " string " << '"' << "positions" << '"' << "" << endl;

    // Teil 7: end remarks
    *Datei << "object " << '"' 
           << "irregular positions irregular connections"
           << '"' << " class field" << 
      endl;
    *Datei << "component " << '"' << "positions" << '"' << " value 1" << endl;
    *Datei << "component " 
           << '"' << "connections" << '"' << " value 2" << endl;
    *Datei << "component " << '"' << "data" << '"' << " value 3" << endl;
    *Datei << "end" << endl; 
  }
}

---