global_discretization_gamma_h_n.cpp

vs.lib is a math library in C++ with a set of linear algebra and integrable / differentiable objects

#include "vs_h.h"
#include "dynamic_array.h"
#include "omega_h.h"
#include "u_h.h"
#include "global_discretization_gamma_h_n.h"

Omega_eh& Gamma_h_i::operator()(int en) {// element selector
	Dynamic_Array_Iterator<Omega_eh> iter(the_omega_eh_array);
   while(iter->element_no() != en && iter.more()) iter++; // searching
   return *(iter.current());
} 

int Gamma_h_i::element_order(int en) const {
	Dynamic_Array_Iterator<Omega_eh> 
			iter((Dynamic_Array<Omega_eh>)the_omega_eh_array);
   while(iter->element_no() != en && iter.more()) iter++; // searching
   return iter.index();
}

C0 Global_Discretization_Gamma_h_i::element_coordinate(int en) {
	int nen = the_gamma_h_i(en).element_node_no(),
   	 nn  = the_gamma_h_i(en).element_node_array()[0],
	    nsd = the_gamma_h_i[nn].no_of_dim();
   C0 ret_val;
   if(nsd > 1) { // build a matrix
		ret_val &= C0(nen, nsd, (double*)NULL);
		for(int i = 0; i < nen; i++) {
			int node_no = the_gamma_h_i(en)[i];
			Node& node = the_gamma_h_i[node_no];
			for(int j = 0; j < nsd; j++)
				ret_val[i][j] = node[j];
		}
   } else { // build a vector
		ret_val &= C0(nen, (double*)NULL);
		for(int i = 0; i < nen; i++) {
			int node_no = the_gamma_h_i(en)[i];
			Node& node = the_gamma_h_i[node_no];
         ret_val[i] = node[0];
		}
   }
	return ret_val;
}

C0 Global_Discretization_Gamma_h_i::element_fixed_variable(int en) {
	int nen = the_gamma_h_i(en).element_node_no(),
   	 nn  = the_gamma_h_i(en).element_node_array()[0],
	    ndf = u_h()[nn].no_of_dim();
	C0 ret_val;
   ret_val &= C0(nen*ndf, (double*)NULL);
   for(int i = 0; i < nen; i++) {
   	int node_no = the_gamma_h_i(en)[i];
      Nodal_Constraint& nc= gh_on_gamma_h()[node_no];
      for(int j = 0; j < ndf; j++)
         if(gh_on_gamma_h().gh_array()[node_no](j)
            == gh_on_Gamma_h::Dirichlet) ret_val[i*ndf+j] = nc[j];
         else ret_val[i*ndf+j] = 0.0;
   }
	return ret_val;
}

C0 Global_Discretization_Gamma_h_i::element_free_variable(int en) {
// residual = - (element_stiffness*element_free_variable)
	int nen = the_gamma_h_i(en).element_node_no(),
   	 nn  = the_gamma_h_i(en).element_node_array()[0],
	    ndf = u_h()[nn].no_of_dim();
	C0 ret_val;
   ret_val &= C0(nen*ndf, (double*)NULL);
   for(int i = 0; i < nen; i++) {
      int node_no = the_gamma_h_i(en)[i];;
      Nodal_Value& var = u_h().u_h_array()[node_no];
      for(int j = 0; j < ndf; j++)
         if(gh_on_gamma_h().gh_array()[node_no](j) == gh_on_Gamma_h::Neumann)
            ret_val[i*ndf+j] = var[j];
   }
	return ret_val;
}

C0 Global_Discretization_Gamma_h_i::element_nodal_force(int en) {
	int nen = the_gamma_h_i(en).element_node_no(),
   	 nn  = the_gamma_h_i(en).element_node_array()[0],
	    ndf = u_h()[nn].no_of_dim();
   C0 ret_val;
   ret_val &= C0(nen*ndf, (double*)NULL);
   for(int i = 0; i < nen; i++) {
      int node_no = the_gamma_h_i(en)[i];
      Nodal_Constraint& nc= gh_on_gamma_h()[node_no];
      for(int j = 0; j < ndf; j++)
         if(gh_on_gamma_h().gh_array()[node_no](j)
            == gh_on_Gamma_h::Neumann) ret_val[i*ndf+j] = nc[j];
         else ret_val[i*ndf+j] = 0.0;
   }
	return ret_val;
}