Flow123d
old_bcd.cc
Go to the documentation of this file.
1 /*
2  * old_bcd.cc
3  *
4  * Created on: Jan 28, 2013
5  * Author: jb
6  */
7 
8 
9 #include "flow/old_bcd.hh"
11 #include "mesh/region.hh"
12 
13 #include "system/tokenizer.hh"
14 #include "boost/lexical_cast.hpp"
15 
17  static OldBcdInput *obcd = new OldBcdInput;
18  return obcd;
19 }
20 /*
21 // set all regions of the given Field<...> @p target
22 template <int spacedim, class Value>
23 void OldBcdInput::set_all( Field<spacedim,Value> &target, const Mesh *mesh) {
24  boost::shared_ptr< FieldElementwise<spacedim, Value> > in_field
25  = boost::make_shared< FieldElementwise<spacedim, Value> >(target.n_comp());
26 
27  target.set_field( mesh->region_db().get_region_set("BOUNDARY"), in_field);
28  target.set_mesh(*mesh);
29 
30 }
31 
32 template <int spacedim, class Value>
33 void OldBcdInput::set_field( Field<spacedim,Value> &target, unsigned int bcd_ele_idx, typename Value::return_type &val) {
34  boost::static_pointer_cast< FieldElementwise<spacedim, Value> >(
35  target[ some_bc_region_ ]
36  )->set_data_row(bcd_ele_idx, val);
37 }
38 */
39 
40 #define DIRICHLET 1
41 #define NEUMANN 2
42 #define NEWTON 3
43 
44 void OldBcdInput::read_flow(const Mesh &mesh, const FilePath &flow_bcd)
45 {
46  using namespace boost;
47 
48  vector< unsigned int *> old_to_new_side_numbering;
49 
50  // in the file the sides are numbered according to opposite nodes as they appear in the MSH file
51  unsigned int sides_0 [1] = {0};
52  old_to_new_side_numbering.push_back( sides_0 );
53  unsigned int sides_1 [2] = {0,1};
54  old_to_new_side_numbering.push_back( sides_1 );
55  unsigned int sides_2 [3] = {0,1,2}; //{0,2,1};
56  old_to_new_side_numbering.push_back( sides_2 );
57  unsigned int sides_3 [4] = {0,1,2,3}; //{3,2,1,0};
58  old_to_new_side_numbering.push_back( sides_3 );
59 
60  // check that all fields has same mesh, reuse it for reader
61  mesh_=&mesh;
62  ASSERT(mesh_ , "Null mesh pointer.\n");
63 
64  /*
65  * - read one flow file, fill fields, make ID list
66  * - read second file, check IDs agains ID list, fill fields
67  */
68 
69  flow_type = std::make_shared< FieldEnum >(1);
70  flow_type->set_mesh(mesh_, true);
71  flow_pressure = std::make_shared< FieldScalar >(1);
72  flow_pressure->set_mesh(mesh_, true);
73  flow_flux = std::make_shared< FieldScalar >(1);
74  flow_flux->set_mesh(mesh_, true);
75  flow_sigma = std::make_shared< FieldScalar >(1);
76  flow_sigma->set_mesh(mesh_, true);
77 
78  Tokenizer tok(flow_bcd);
79  try {
80  double scalar, flux, sigma;
81  unsigned int id;
82 
83  xprintf(MsgLog, "Reading old BCD file for flow: %s ...", tok.f_name().c_str());
84  tok.skip_to("$BoundaryConditions");
85  tok.next_line(false);
86  unsigned int n_boundaries = lexical_cast<unsigned int>(*tok); ++tok;
87 
88  for(unsigned int i_bcd=0; i_bcd < n_boundaries; i_bcd++) {
89  tok.next_line();
90 
91  id = lexical_cast<unsigned int>(*tok); ++tok;
92 
93  unsigned int type = lexical_cast<unsigned int>(*tok); ++tok;
94 
95  switch( type ) {
96  case DIRICHLET:
97  scalar = lexical_cast<double>(*tok); ++tok;
98  flux = 0.0;
99  sigma = 0.0;
100  break;
101  case NEUMANN:
102  flux = lexical_cast<double>(*tok); ++tok;
103  sigma = 0.0;
104  scalar = 0.0;
105  break;
106  case NEWTON:
107  scalar = lexical_cast<double>(*tok); ++tok;
108  sigma = lexical_cast<double>(*tok); ++tok;
109  flux = 0.0;
110  break;
111  default :
112  xprintf(UsrErr,"Unknown type of boundary condition - cond # %d, type %c\n", id, type );
113  break;
114  }
115 
116  unsigned int where = lexical_cast<unsigned int>(*tok); ++tok;
117 
118  unsigned int eid, sid, bc_ele_idx, our_sid;
119  Element * ele;
120  Boundary * bcd;
121 
122  switch( where ) {
123  case 2: // SIDE_EL - BC given by element and its local side number
124  eid = lexical_cast<unsigned int>(*tok); ++tok;
125  sid = lexical_cast<unsigned int>(*tok); ++tok;
126 
127  // find and set the side
128  // const cast can be removed when get rid of FullIterators and whole sys_vector stuff
129  // and have correct constantness for mesh classes
130  ele = const_cast<Element *>(mesh_->element.find_id( eid ));
131  if( sid < 0 || sid >= ele->n_sides() )
132  xprintf(UsrErr,"Boundary %d has incorrect reference to side %d\n", id, sid );
133  our_sid=old_to_new_side_numbering[ele->dim()][sid];
134  bcd = ele->side(our_sid) -> cond();
135  if (bcd) {
136  bc_ele_idx = mesh_->bc_elements.index( ele->side(our_sid) -> cond()->element() );
137  id_2_bcd_[id]= bc_ele_idx;
138  if ( ! some_bc_region_.is_valid() ) some_bc_region_ = ele->side(our_sid) -> cond()->element()->region();
139 
140  flow_type->set_data_row( bc_ele_idx, type);
141  flow_pressure->set_data_row(bc_ele_idx, scalar);
142  flow_flux->set_data_row( bc_ele_idx, flux);
143  flow_sigma->set_data_row( bc_ele_idx, sigma);
144  } else {
145  xprintf(Warn, "IGNORING boundary condition %d for non-boundary side %d of element ID: %d\n", id, sid, eid);
146  }
147  break;
148  case 3: // SIDE_E - BC given only by element, apply to all its sides
149 
150  xprintf(UsrErr, "Element only BCD are not supported.\n");
151  /*
152  eid = atoi( xstrtok( NULL) );
153 
154  // find and set all exterior sides, possibly add more boundaries
155  ele = mesh->element.find_id( eid );
156  n_exterior=0;
157  FOR_ELEMENT_SIDES(ele, li) {
158  sde = ele->side( li );
159  if ( bcd=sde->cond() ) {
160 
161  if (n_exterior > 0) {
162  xprintf(UsrErr, "Implicitly setting BC %d on more then one exterior sides of the element %d.\n", bcd_id, eid);
163  //BoundaryFullIter new_bcd = mesh->boundary.add_item();
164  // *new_bcd = *bcd;
165  //bcd=new_bcd;
166  }
167  bcd->type = type;
168  bcd->flux = flux;
169  bcd->scalar = scalar;
170  bcd->sigma = sigma;
171  n_exterior++;
172  }
173  }
174  */
175  break;
176  default:
177  xprintf(UsrErr,"Unknown entity for boundary condition - cond # %d, ent. %c\n", id, where );
178  break;
179  }
180  unsigned int n_tags = lexical_cast<unsigned int>(*tok); ++tok;
181  while (n_tags>0) ++tok, --n_tags; // skip remaining tags
182 
183 
184 
185  // There was possibility to set group IDs for boundary faces (like regions for boundary elements)
186  // It is probably not used so we do not implement it as it is DEPRECATED
187  /*
188  n_tags = atoi( xstrtok( NULL) );
189  if( n_tags > 0 ) {
190  int group_id = atoi( xstrtok( NULL) );
191  flow::VectorId<int>::FullIter group_iter( mesh->bcd_group_id.find_id(group_id) );
192 
193  if ( group_iter == mesh->bcd_group_id.end() ) {
194  // not found -> create new group
195  group_iter = mesh->bcd_group_id.add_item(group_id);
196  }
197  bcd->group = group_iter.index(); // in fact we do not use integres stored in the vector, but we use index
198  }
199  */
200  }
201  xprintf(MsgLog, "DONE\n");
202  } catch (bad_lexical_cast &) {
203  xprintf(UsrErr, "Wrong format of number, %s.\n", tok.position_msg().c_str());
204  } // flow bcd reader
205 }
206 
207 
208 
209 void OldBcdInput::read_transport(unsigned int n_substances, const FilePath &transport_bcd)
210 {
211  using namespace boost;
212 
213  ASSERT(mesh_ , "Null mesh pointer.\n");
214  trans_conc = std::make_shared< FieldVector >( n_substances );
215  trans_conc->set_mesh(mesh_, true);
216 
217  FieldValue<3>::Vector::return_type ele_value(n_substances);
218 
219  Tokenizer tok(transport_bcd);
220  try {
221  unsigned int bcd_id, boundary_id, bc_ele_idx;
222 
223  xprintf(MsgLog, "Reading old BCD file for transport: %s ...", tok.f_name().c_str());
224  if (tok.skip_to("$Transport_BCDFormat")) tok.next_line(false);
225  tok.skip_to("$Transport_BCD");
226  tok.next_line(false);
227  unsigned int n_bcd = lexical_cast<unsigned int>(*tok); ++tok;
228  for (unsigned int i_bcd = 0; i_bcd < n_bcd; i_bcd++) {
229  tok.next_line();
230  bcd_id = lexical_cast<unsigned int>(*tok); ++tok;
231  boundary_id = lexical_cast<unsigned int>(*tok); ++tok;
232 
234  if (it == id_2_bcd_.end())
235  xprintf(UsrErr,"Wrong boundary index %d for bcd id %d in transport bcd file!", boundary_id, bcd_id);
236  bc_ele_idx = it->second;
237 
238  for (unsigned int sbi = 0; sbi < n_substances; sbi++) {
239  ele_value[sbi] = lexical_cast<double>(*tok); ++tok;
240  }
241 
242  trans_conc->set_data_row(bc_ele_idx, ele_value);
243 
244  }
245 
246  xprintf(MsgLog, "DONE\n");
247 
248  } catch (bad_lexical_cast &) {
249  xprintf(UsrErr, "Wrong format of number, %s.\n", tok.position_msg().c_str());
250  } // flow bcd reader
251 }