LCOV - code coverage report
Current view: top level - ugbase/lib_grid/file_io - file_io_grdecl.cpp (source / functions) Coverage Total Hit
Test: coverage.info Lines: 0.0 % 274 0
Test Date: 2025-09-21 23:31:46 Functions: 0.0 % 8 0

            Line data    Source code
       1              : /*
       2              :  * Copyright (c) 2009-2015:  G-CSC, Goethe University Frankfurt
       3              :  * Author: Shuai Lu
       4              :  * 
       5              :  * This file is part of UG4.
       6              :  * 
       7              :  * UG4 is free software: you can redistribute it and/or modify it under the
       8              :  * terms of the GNU Lesser General Public License version 3 (as published by the
       9              :  * Free Software Foundation) with the following additional attribution
      10              :  * requirements (according to LGPL/GPL v3 §7):
      11              :  * 
      12              :  * (1) The following notice must be displayed in the Appropriate Legal Notices
      13              :  * of covered and combined works: "Based on UG4 (www.ug4.org/license)".
      14              :  * 
      15              :  * (2) The following notice must be displayed at a prominent place in the
      16              :  * terminal output of covered works: "Based on UG4 (www.ug4.org/license)".
      17              :  * 
      18              :  * (3) The following bibliography is recommended for citation and must be
      19              :  * preserved in all covered files:
      20              :  * "Reiter, S., Vogel, A., Heppner, I., Rupp, M., and Wittum, G. A massively
      21              :  *   parallel geometric multigrid solver on hierarchically distributed grids.
      22              :  *   Computing and visualization in science 16, 4 (2013), 151-164"
      23              :  * "Vogel, A., Reiter, S., Rupp, M., Nägel, A., and Wittum, G. UG4 -- a novel
      24              :  *   flexible software system for simulating pde based models on high performance
      25              :  *   computers. Computing and visualization in science 16, 4 (2013), 165-179"
      26              :  * 
      27              :  * This program is distributed in the hope that it will be useful,
      28              :  * but WITHOUT ANY WARRANTY; without even the implied warranty of
      29              :  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
      30              :  * GNU Lesser General Public License for more details.
      31              :  */
      32              : #include <iostream>
      33              : #include <fstream>
      34              : #include <sstream>
      35              : #include <string>
      36              : #include <vector>
      37              : 
      38              : #include "../lg_base.h"
      39              : #include "common/util/string_util.h"
      40              : #include "lib_grid/global_attachments.h"
      41              : 
      42              : #include "file_io_grdecl.h"
      43              : 
      44              : using namespace std;
      45              : 
      46              : namespace ug
      47              : {
      48              : 
      49              :         struct xy
      50              :         {
      51              :            double x, y;
      52              :         };
      53              :         
      54              :         struct xyz
      55              :         {
      56              :            double x, y, z;
      57              :         };
      58              : 
      59              :         struct ab
      60              :         {
      61              :                 size_t a1, a2, a3, a4, b1, b2, b3, b4;
      62              :         };
      63              : 
      64              : 
      65            0 : void GetDim(string str, vector<size_t>& dim)
      66              : {
      67            0 :     string word = "";
      68            0 :     for (auto x : str) 
      69              :     {
      70            0 :         if (x == ' ')
      71              :         {
      72            0 :                         if (word!="" && dim.size()<3)
      73              :                         {
      74              :                                 size_t a;
      75            0 :                                 stringstream ss;
      76              :                                 ss << word;
      77              :                                 ss >> a;
      78            0 :                                 dim.push_back(a);
      79            0 :                         }  
      80              :             word = "";
      81              :         }
      82              :         else {
      83            0 :             word = word + x;
      84              :         }
      85              :     }
      86            0 :         if (word!="" && dim.size()<3)
      87              :                 {
      88              :                         size_t a;
      89            0 :                         stringstream ss;
      90              :                         ss << word;
      91              :                         ss >> a;
      92            0 :                         dim.push_back(a);
      93            0 :                 }
      94            0 : }
      95              : 
      96            0 : void GetCoord(string str, vector<double>& coord)
      97              : {
      98            0 :     string word = "";
      99            0 :     for (auto x : str) 
     100              :     {
     101            0 :         if (x == ' ')
     102              :         {
     103            0 :                         if (word!="")
     104              :                         {
     105            0 :                                 string snum = "";
     106            0 :                                 size_t num=1;
     107            0 :                                 for (auto y : word)
     108              :                                 {
     109            0 :                                         if (y=='*')
     110              :                                         {
     111            0 :                                                 stringstream ss;
     112              :                                                 ss << snum;
     113              :                                                 ss >> num;
     114              :                                                 snum="";
     115            0 :                                         }
     116              :                                         else
     117            0 :                                                 snum=snum+y;
     118              :                                 }
     119              :                                 
     120              :                                 double value;
     121            0 :                                 stringstream ss;
     122              :                                 ss << snum;
     123              :                                 ss >> value;
     124            0 :                                 for (size_t i=0; i<num; i++)
     125              :                                 {
     126            0 :                                         coord.push_back(value);
     127              :                                 }
     128            0 :                         }  
     129              :             word = "";
     130              :         }
     131              :         else {
     132            0 :             word = word + x;
     133              :         }
     134              :     }
     135            0 :         if (word!="" && word.compare(0,1,"/")!=0 && word!="\r")
     136              :         {
     137            0 :                 string snum = "";
     138            0 :                 size_t num=1;
     139            0 :                 for (auto y : word)
     140              :                 {
     141            0 :                         if (y=='*')
     142              :                         {
     143            0 :                                 stringstream ss;
     144              :                                 ss << snum;
     145              :                                 ss >> num;
     146              :                                 snum="";
     147            0 :                         }
     148              :                         else
     149            0 :                                 snum=snum+y;
     150              :                 }
     151              :                                 
     152              :                 double value;
     153            0 :                 stringstream ss;
     154              :                 ss << snum;
     155              :                 ss >> value;
     156            0 :                 for (size_t i=0; i<num; i++)
     157              :                 {
     158            0 :                         coord.push_back(value);
     159              :                 }
     160            0 :         }
     161            0 : }
     162              :         
     163            0 : void GetZcorn(string str, vector<double>& zcorn)
     164              : {
     165            0 :     string word = "";
     166            0 :     for (auto x : str) 
     167              :     {
     168            0 :         if (x == ' ')
     169              :         {
     170            0 :                         if (word!="")
     171              :                         {
     172            0 :                                 string snum = "";
     173            0 :                                 size_t num=1;
     174            0 :                                 for (auto y : word)
     175              :                                 {
     176            0 :                                         if (y=='*')
     177              :                                         {
     178            0 :                                                 stringstream ss;
     179              :                                                 ss << snum;
     180              :                                                 ss >> num;
     181              :                                                 snum="";
     182            0 :                                         }
     183              :                                         else
     184            0 :                                                 snum=snum+y;
     185              :                                 }
     186              :                                 
     187              :                                 double value;
     188            0 :                                 stringstream ss;
     189              :                                 ss << snum;
     190              :                                 ss >> value;
     191            0 :                                 for (size_t i=0; i<num; i++)
     192              :                                 {
     193            0 :                                         zcorn.push_back(value);
     194              :                                 }
     195            0 :                         }  
     196              :             word = "";
     197              :         }
     198              :         else {
     199            0 :             word = word + x;
     200              :         }
     201              :     }
     202            0 :         if (word!="" && word.compare(0,1,"/")!=0)
     203              :         {
     204            0 :                 string snum = "";
     205            0 :                 size_t num=1;
     206            0 :                 for (auto y : word)
     207              :                 {
     208            0 :                         if (y=='*')
     209              :                         {
     210            0 :                                 stringstream ss;
     211              :                                 ss << snum;
     212              :                                 ss >> num;
     213              :                                 snum="";
     214            0 :                         }
     215              :                         else
     216            0 :                                 snum=snum+y;
     217              :                 }
     218              :                                 
     219              :                 double value;
     220            0 :                 stringstream ss;
     221              :                 ss << snum;
     222              :                 ss >> value;
     223            0 :                 for (size_t i=0; i<num; i++)
     224              :                 {
     225            0 :                         zcorn.push_back(value);
     226              :                 }
     227            0 :         }
     228            0 : }
     229              : 
     230            0 : void GetAct(string str, vector<size_t>& act)
     231              : {
     232            0 :     string word = "";
     233            0 :     for (auto x : str) 
     234              :     {
     235            0 :         if (x == ' ')
     236              :         {
     237            0 :                         if (word!="")
     238              :                         {
     239            0 :                                 string snum = "";
     240            0 :                                 size_t num=1;
     241            0 :                                 for (auto y : word)
     242              :                                 {
     243            0 :                                         if (y=='*')
     244              :                                         {
     245            0 :                                                 stringstream ss;
     246              :                                                 ss << snum;
     247              :                                                 ss >> num;
     248              :                                                 snum="";
     249            0 :                                         }
     250              :                                         else
     251            0 :                                                 snum=snum+y;
     252              :                                 }
     253              :                                 
     254              :                                 size_t value;
     255            0 :                                 stringstream ss;
     256              :                                 ss << snum;
     257              :                                 ss >> value;
     258            0 :                                 for (size_t i=0; i<num; i++)
     259              :                                 {
     260            0 :                                         act.push_back(value);
     261              :                                 }
     262            0 :                         }  
     263              :             word = "";
     264              :         }
     265              :         else {
     266            0 :             word = word + x;
     267              :         }
     268              :     }
     269            0 :         if (word!="" && word.compare(0,1,"/")!=0)
     270              :         {
     271            0 :                 string snum = "";
     272            0 :                 size_t num=1;
     273            0 :                 for (auto y : word)
     274              :                 {
     275            0 :                         if (y=='*')
     276              :                         {
     277            0 :                                 stringstream ss;
     278              :                                 ss << snum;
     279              :                                 ss >> num;
     280              :                                 snum="";
     281            0 :                         }
     282              :                         else
     283            0 :                                 snum=snum+y;
     284              :                 }
     285              :                                 
     286              :                 size_t value;
     287            0 :                 stringstream ss;
     288              :                 ss << snum;
     289              :                 ss >> value;
     290            0 :                 for (size_t i=0; i<num; i++)
     291              :                 {
     292            0 :                         act.push_back(value);
     293              :                 }
     294            0 :         }
     295            0 : }
     296              : 
     297            0 : void GetProperty(string str, vector<double>& prop)
     298              : {
     299            0 :     string word = "";
     300            0 :     for (auto x : str) 
     301              :     {
     302            0 :         if (x == ' ')
     303              :         {
     304            0 :                         if (word!="")
     305              :                         {
     306            0 :                                 string snum = "";
     307            0 :                                 size_t num=1;
     308            0 :                                 for (auto y : word)
     309              :                                 {
     310            0 :                                         if (y=='*')
     311              :                                         {
     312            0 :                                                 stringstream ss;
     313              :                                                 ss << snum;
     314              :                                                 ss >> num;
     315              :                                                 snum="";
     316            0 :                                         }
     317              :                                         else
     318            0 :                                                 snum=snum+y;
     319              :                                 }
     320              :                                 
     321              :                                 double value;
     322            0 :                                 stringstream ss;
     323              :                                 ss << snum;
     324              :                                 ss >> value;
     325            0 :                                 for (size_t i=0; i<num; i++)
     326              :                                 {
     327            0 :                                         prop.push_back(value);
     328              :                                 }
     329            0 :                         }  
     330              :             word = "";
     331              :         }
     332              :         else {
     333            0 :             word = word + x;
     334              :         }
     335              :     }
     336            0 :         if (word!="" && word.compare(0,1,"/")!=0)
     337              :         {
     338            0 :                 string snum = "";
     339            0 :                 size_t num=1;
     340            0 :                 for (auto y : word)
     341              :                 {
     342            0 :                         if (y=='*')
     343              :                         {
     344            0 :                                 stringstream ss;
     345              :                                 ss << snum;
     346              :                                 ss >> num;
     347              :                                 snum="";
     348            0 :                         }
     349              :                         else
     350            0 :                                 snum=snum+y;
     351              :                 }
     352              :                                 
     353              :                 double value;
     354            0 :                 stringstream ss;
     355              :                 ss << snum;
     356              :                 ss >> value;
     357            0 :                 for (size_t i=0; i<num; i++)
     358              :                 {
     359            0 :                         prop.push_back(value);
     360              :                 }
     361            0 :         }
     362            0 : }
     363              : 
     364              : 
     365            0 : bool AttachAct(Grid& grid, const char* filename, string name)
     366              : {
     367              :         vector<size_t> actnum;
     368              :         bool d = false;
     369              :         string buf;
     370              :         //string name2 = name +' ';
     371            0 :         ifstream ifs(filename);
     372            0 :         if(!ifs)
     373              :                 return false;
     374              : 
     375            0 :         while (getline(ifs, buf))
     376              :         {
     377            0 :                 if (d==true)
     378              :                         {
     379            0 :                                 if (buf.compare(2,1,"/")==0)
     380              :                                         d=false;
     381            0 :                                 else if (buf.compare(buf.size()-2,1,"/")==0)
     382              :                                 {
     383            0 :                                         GetAct(buf, actnum);
     384              :                                         d=false;
     385              :                                 }
     386              :                                 else
     387            0 :                                         GetAct(buf, actnum);
     388              :                         }
     389            0 :                 if (buf.compare(0,name.length(),name)==0)
     390              :                         d=true;
     391              :         }
     392              :         
     393            0 :         if (!GlobalAttachments::is_declared(name))
     394              :         {
     395            0 :                 try {GlobalAttachments::declare_attachment<ANumber>(name);}
     396            0 :                 catch (ug::UGError& err)
     397              :                 {
     398              :                         UG_LOGN(err.get_msg());
     399              :                         return false;
     400            0 :                 }
     401              :         }
     402            0 :         ANumber aVols = GlobalAttachments::attachment<ANumber>(name);
     403            0 :         if(!grid.has_volume_attachment(aVols))
     404              :                 grid.attach_to_volumes(aVols);
     405              : 
     406              :         Grid::VolumeAttachmentAccessor<ANumber> aaVols(grid, aVols);
     407              :         
     408              :         size_t i = 0;
     409            0 :         for(VolumeIterator iter = grid.volumes_begin(); iter != grid.volumes_end(); ++iter, i++)
     410              :                 {
     411            0 :                         aaVols[*iter]=actnum[i];
     412              : //                      ofs<<i<<' '<<aaVols[*iter]<<endl;
     413              :                 }
     414              :         return true;
     415            0 : }
     416              : 
     417              : 
     418            0 : bool AttachProperty(Grid& grid, const char* filename, string name)
     419              : {
     420              :         vector<double> prop;
     421              :         bool d = false;
     422              :         bool e = false;
     423              :         string buf;
     424            0 :         string name2 = name +' ';
     425            0 :         ifstream ifs(filename);
     426            0 :         if(!ifs)
     427              :                 return false;
     428              : 
     429            0 :         while (getline(ifs, buf))
     430              :         {
     431            0 :                 if (d==true)
     432              :                         {
     433            0 :                                 if (buf.compare(2,1,"/")==0)
     434              :                                         d=false;
     435            0 :                                 else if (buf.compare(buf.size()-2,1,"/")==0)
     436              :                                 {
     437            0 :                                         GetProperty(buf, prop);
     438              :                                         d=false;
     439              :                                 }
     440              :                                 else
     441            0 :                                         GetProperty(buf, prop);
     442              :                         }
     443            0 :                 if ((buf.compare(0,11,"-- Property")==0)&&(e==true))
     444              :                         d=true;
     445            0 :                 if (buf.compare(0,name2.length(),name2)==0)
     446              :                         e=true;
     447              :         }
     448              :         
     449            0 :         if (!GlobalAttachments::is_declared(name))
     450              :         {
     451            0 :                 try {GlobalAttachments::declare_attachment<ANumber>(name);}
     452            0 :                 catch (ug::UGError& err)
     453              :                 {
     454              :                         UG_LOGN(err.get_msg());
     455              :                         return false;
     456            0 :                 }
     457              :         }
     458            0 :         ANumber aVols = GlobalAttachments::attachment<ANumber>(name);
     459            0 :         if(!grid.has_volume_attachment(aVols))
     460              :                 grid.attach_to_volumes(aVols);
     461              : 
     462              :         Grid::VolumeAttachmentAccessor<ANumber> aaVols(grid, aVols);
     463              :         
     464              : //      ofstream ofs;
     465              : //      ofs.open ("Volumes.txt");
     466              :         
     467              :         size_t i = 0;
     468            0 :         for(VolumeIterator iter = grid.volumes_begin(); iter != grid.volumes_end(); ++iter, i++)
     469              :                 {
     470            0 :                         aaVols[*iter]=prop[i];
     471              : //                      ofs<<i<<' '<<aaVols[*iter]<<endl;
     472              :                 }
     473              :         return true;
     474            0 : }
     475              : 
     476            0 : bool LoadGridFromGRDECL(Grid& grid, const char* filename, AVector3& aPos)
     477              : {
     478              :         
     479              :         string buf;
     480              :         vector<size_t> dim;
     481              :         vector<double> coord;
     482              :         vector<double> zcorn;
     483              :         vector<size_t> act;
     484              :         vector<string> prop;
     485              :         vector<string> prop_name;
     486              :         
     487              :         bool a=false;
     488              :         bool b=false;
     489              :         bool c=false;
     490              :         bool d=false;
     491              :         
     492            0 :         ifstream ifs(filename);
     493            0 :         if(!ifs)
     494              :                 return false;
     495              :         
     496            0 :         while (getline(ifs, buf))
     497              :                 {
     498              :                         //Get the dimension of the cells
     499            0 :                         if (a==true)
     500              :                         {
     501            0 :                                 GetDim(buf, dim);
     502              :                                 a=false;
     503              :                         }
     504            0 :                         if (buf.compare(0,8,"SPECGRID")==0)
     505              :                                 a=true;
     506              :                         
     507              :                         //Get the coordinate list of the top and bot surfaces
     508            0 :                         if (b==true)
     509              :                                 {
     510            0 :                                         if (buf.compare(2,1,"/")==0)
     511              :                                                 b=false;
     512            0 :                                         else if (buf.compare(buf.size()-2,1,"/")==0)
     513              :                                         {
     514            0 :                                                 GetCoord(buf, coord);
     515              :                                                 b=false;
     516              :                                         }
     517              :                                         else
     518            0 :                                                 GetCoord(buf, coord);
     519              :                                 }
     520            0 :                         if (buf.compare(0,6,"COORD ")==0 || buf.compare(0,7,"COORD\r")==0)
     521              :                                 b=true;
     522              :                                 
     523              :                         //Get the depth list
     524            0 :                         if (c==true)
     525              :                                 {
     526            0 :                                         if (buf.compare(2,1,"/")==0)
     527              :                                                 c=false;
     528            0 :                                         else if (buf.compare(buf.size()-2,1,"/")==0)
     529              :                                         {
     530            0 :                                                 GetCoord(buf, zcorn);
     531              :                                                 c=false;
     532              :                                         }
     533              :                                         else
     534            0 :                                                 GetCoord(buf, zcorn);
     535              :                                 }
     536            0 :                         if (buf.compare(0,5,"ZCORN")==0)
     537              :                                 c=true;
     538              :                                 
     539              :                         //Get the ACTNUM list
     540            0 :                         if (d==true)
     541              :                                 {
     542            0 :                                         if (buf.compare(2,1,"/")==0)
     543              :                                                 d=false;
     544            0 :                                         else if (buf.compare(buf.size()-2,1,"/")==0)
     545              :                                         {
     546            0 :                                                 GetAct(buf, act);
     547              :                                                 d=false;
     548              :                                         }
     549              :                                         else
     550            0 :                                                 GetAct(buf, act);
     551              :                                 }
     552            0 :                         if (buf.compare(0,6,"ACTNUM")==0)
     553              :                                 d=true;
     554              :         
     555              :                         //Get property name list
     556            0 :                         prop.push_back(buf);
     557            0 :                         if (buf.compare(0,11,"-- Property")==0)
     558              :                         {
     559              :                                 prop.pop_back();
     560            0 :                                 string name = "";
     561            0 :                                 for (auto x : prop.back()) 
     562              :                                 {
     563            0 :                                         if (x != ' ')
     564            0 :                                                 name = name + x;
     565              :                                         else
     566              :                                                 break;
     567              :                                 }
     568            0 :                                 prop_name.push_back(name);
     569              :                         }
     570              :                         
     571              :                 }
     572              : 
     573              :                 
     574              :         vector<xy> top;
     575              :         vector<xy> bot;
     576              :         // loop over pillar
     577            0 :         for (size_t j=0; j<dim[1]+1; j++)
     578              :         {
     579            0 :                 for (size_t i=0; i<dim[0]+1; i++)
     580              :                 {
     581              :                         //build pillar coords
     582            0 :                         top.push_back({coord[(j*(dim[0]+1)+i)*6], coord[(j*(dim[0]+1)+i)*6+1]});
     583            0 :                         bot.push_back({coord[(j*(dim[0]+1)+i)*6+3], coord[(j*(dim[0]+1)+i)*6+4]});
     584              :                         //check pillar is vertical
     585            0 :                         if ( coord[(j*(dim[0]+1)+i)*6] != coord[(j*(dim[0]+1)+i)*6+3]|| coord[(j*(dim[0]+1)+i)*6+1] != coord[(j*(dim[0]+1)+i)*6+4] )
     586            0 :                                 UG_THROW ("LoadGridFromGRDECL: pillar is not vertical");
     587              :                 }
     588              :         }
     589              :         
     590            0 :         vector<vector<double> > zcorn_ij((dim[0]+1)*(dim[1]+1), vector<double>(1,0));
     591              :         
     592            0 :         for (size_t k=0; k<dim[2]; k++)
     593              :         {
     594            0 :                 for (size_t j=0; j<dim[1]; j++)
     595              :                 {
     596            0 :                         for (size_t i=0; i<dim[0]; i++)
     597              :                         {
     598            0 :                                 size_t ij=j*(dim[0]+1)+i;
     599            0 :                                 size_t ijk4=k*dim[1]*dim[0]*8+j*dim[0]*4+i*2;
     600              :                                 
     601            0 :                                 zcorn_ij[ij].push_back(zcorn[ijk4]);
     602            0 :                                 zcorn_ij[ij+1].push_back(zcorn[ijk4+1]);
     603            0 :                                 zcorn_ij[ij+dim[0]+2].push_back(zcorn[ijk4+dim[0]*2+1]);
     604            0 :                                 zcorn_ij[ij+dim[0]+1].push_back(zcorn[ijk4+dim[0]*2]);
     605              :                                 
     606            0 :                                 zcorn_ij[ij].push_back(zcorn[ijk4+dim[1]*dim[0]*4]);
     607            0 :                                 zcorn_ij[ij+1].push_back(zcorn[ijk4+dim[1]*dim[0]*4+1]);
     608            0 :                                 zcorn_ij[ij+dim[0]+2].push_back(zcorn[ijk4+dim[1]*dim[0]*4+dim[0]*2+1]);
     609            0 :                                 zcorn_ij[ij+dim[0]+1].push_back(zcorn[ijk4+dim[1]*dim[0]*4+dim[0]*2]);                  
     610              :                         }
     611              :                 }
     612              :         }
     613              :         
     614            0 :         vector<vector<size_t> > Index_zcorn((dim[0]+1)*(dim[1]+1), vector<size_t>(1,0));
     615            0 :         vector<vector<size_t> > Index_zcorn_Global((dim[0]+1)*(dim[1]+1), vector<size_t>(1,0));
     616            0 :         vector<vector<double> > New_zcorn((dim[0]+1)*(dim[1]+1), vector<double>(1,0));
     617            0 :         size_t numVrts=0;
     618              :         
     619              :         vector<xyz> coord_list;
     620              : 
     621            0 :         for (size_t j=0; j<dim[1]+1; j++)
     622              :         {
     623            0 :                 for (size_t i=0; i<dim[0]+1; i++)
     624              :                 {       
     625            0 :                         size_t ij=j*(dim[0]+1)+i;
     626              :                         
     627              :                         //remove redundant coords
     628            0 :                         New_zcorn[ij].push_back(zcorn_ij[ij][1]);
     629            0 :                         for (size_t t=2; t<zcorn_ij[ij].size(); t++)
     630              :                         {
     631              :                                 bool New_coord=true;
     632            0 :                                 for (size_t s =1; s<New_zcorn[ij].size(); s++)
     633              :                                 {
     634            0 :                                         if (zcorn_ij[ij][t]==New_zcorn[ij][s])
     635              :                                         {
     636              :                                                 New_coord = false;
     637            0 :                                                 Index_zcorn[ij].push_back(s-1);
     638              :                                                 break;
     639              :                                         }
     640              :                                 }
     641              :                                 
     642              :                                 if (New_coord==true)
     643              :                                 {
     644            0 :                                         Index_zcorn[ij].push_back(New_zcorn[ij].size()-1);
     645            0 :                                         New_zcorn[ij].push_back(zcorn_ij[ij][t]);
     646              :                                 }
     647              :                         }
     648              :                         
     649            0 :                         for (size_t k=1; k<New_zcorn[ij].size(); k++)
     650              :                         {
     651            0 :                                 coord_list.push_back({top[ij].x, top[ij].y, -New_zcorn[ij][k]});
     652            0 :                                 Index_zcorn_Global[ij].push_back(numVrts);
     653            0 :                                 numVrts++;
     654              :                         }
     655              :                         
     656              :                         reverse(New_zcorn[ij].begin(), New_zcorn[ij].end());
     657              :                         New_zcorn[ij].pop_back();
     658              :                         reverse(zcorn_ij[ij].begin(), zcorn_ij[ij].end());
     659              :                         zcorn_ij[ij].pop_back();
     660              :                         reverse(Index_zcorn[ij].begin(), Index_zcorn[ij].end());
     661              :                         
     662              :                         reverse(Index_zcorn_Global[ij].begin(), Index_zcorn_Global[ij].end());
     663              :                         Index_zcorn_Global[ij].pop_back();
     664              :                         reverse(Index_zcorn_Global[ij].begin(), Index_zcorn_Global[ij].end());
     665              :                 }
     666              :                 
     667              :         }
     668              :         
     669              :         vector<ab> ele_list;
     670              :         
     671            0 :         for (size_t k=0; k<dim[2]; k++)
     672              :         {
     673            0 :                 for (size_t j=0; j<dim[1]; j++)
     674              :                 {
     675            0 :                         for (size_t i=0; i<dim[0]; i++)
     676              :                         {
     677            0 :                                 size_t ij=j*(dim[0]+1)+i;
     678              :                                 
     679            0 :                                 size_t a1=Index_zcorn_Global[ij][Index_zcorn[ij].back()];
     680              :                                 Index_zcorn[ij].pop_back();
     681            0 :                                 size_t a2=Index_zcorn_Global[ij+1][Index_zcorn[ij+1].back()];
     682              :                                 Index_zcorn[ij+1].pop_back();
     683            0 :                                 size_t a3=Index_zcorn_Global[ij+dim[0]+2][Index_zcorn[ij+dim[0]+2].back()];
     684              :                                 Index_zcorn[ij+dim[0]+2].pop_back();
     685            0 :                                 size_t a4=Index_zcorn_Global[ij+dim[0]+1][Index_zcorn[ij+dim[0]+1].back()];
     686              :                                 Index_zcorn[ij+dim[0]+1].pop_back();
     687              :                                 
     688            0 :                                 size_t b1=Index_zcorn_Global[ij][Index_zcorn[ij].back()];
     689              :                                 Index_zcorn[ij].pop_back();
     690            0 :                                 size_t b2=Index_zcorn_Global[ij+1][Index_zcorn[ij+1].back()];
     691              :                                 Index_zcorn[ij+1].pop_back();
     692            0 :                                 size_t b3=Index_zcorn_Global[ij+dim[0]+2][Index_zcorn[ij+dim[0]+2].back()];
     693              :                                 Index_zcorn[ij+dim[0]+2].pop_back();
     694            0 :                                 size_t b4=Index_zcorn_Global[ij+dim[0]+1][Index_zcorn[ij+dim[0]+1].back()];
     695              :                                 Index_zcorn[ij+dim[0]+1].pop_back();
     696              :                                 
     697            0 :                                 ele_list.push_back({a1, a2, a3, a4, b1, b2, b3, b4});
     698              :                         }
     699              :                 }
     700              :         }
     701              :                 
     702              :         size_t numElems = ele_list.size();
     703              : 
     704              : //      create points
     705              : //      store pointers to the vertices on the fly in a vector.
     706              :         vector<Vertex*>   vVrts;
     707              :         vector<size_t> vVrtIds;
     708            0 :         vVrts.resize(numVrts); vVrtIds.resize(numVrts);
     709              : 
     710            0 :         for(size_t i = 0; i < numVrts; ++i)
     711            0 :                 vVrts[i] = *grid.create<RegularVertex>();
     712              : 
     713            0 :         if(!grid.has_vertex_attachment(aPos))
     714            0 :                 grid.attach_to_vertices(aPos);
     715              :         Grid::VertexAttachmentAccessor<APosition> aaPos(grid, aPos);
     716              : 
     717              : //      read the points
     718              :         {
     719              :                 size_t i = 0;
     720            0 :                 for(VertexIterator iter = grid.vertices_begin(); iter != grid.vertices_end(); ++iter, ++i)
     721              :                 {
     722            0 :                         vVrtIds[i]=i;
     723            0 :                         aaPos[*iter].x()=coord_list[i].x;
     724            0 :                         aaPos[*iter].y()=coord_list[i].y;
     725            0 :                         aaPos[*iter].z()=coord_list[i].z;
     726              :                 }
     727              : 
     728              :         }
     729              : 
     730              : //      read the hexahedrons
     731              :         {
     732            0 :                 for(size_t i = 0; i < numElems; ++i)
     733              :                 {
     734              :                         grid.create<Hexahedron>
     735            0 :                                 (HexahedronDescriptor
     736            0 :                                         (vVrts[ele_list[i].a1], vVrts[ele_list[i].a2], vVrts[ele_list[i].a3], vVrts[ele_list[i].a4], vVrts[ele_list[i].b1], vVrts[ele_list[i].b2], vVrts[ele_list[i].b3], vVrts[ele_list[i].b4]));
     737              :                 }
     738              :         }
     739              : 
     740            0 :         AttachAct(grid, filename, "ACTNUM");  
     741              : 
     742              : //      read the volumes(properties)
     743              :         
     744            0 :         for (size_t i = 0; i<prop_name.size(); i++)
     745              :                 {
     746            0 :                         AttachProperty(grid, filename, prop_name[i]);
     747              :                 }
     748              : 
     749              : 
     750              :         return true;
     751            0 : }
     752              : 
     753              : }//     end of namespace
        

Generated by: LCOV version 2.0-1