wztensorgrid.cxx

有限元学习研究用源代码(老外的),供科研人员参考

#include "wztensorgrid.hxx"

wzFloat wzTensorGrid::epsilon = 1.e-7;

wzTensorGrid::wzTensorGrid(wzIndex dim, wzFloat eps)
{
  int i;
  for(i=0;i<wzPointDim;i++){
    x[i].reinitialize(n[i]);
    r[i] = 0;
  }
  checked = 0;
}
void wzTensorGrid::check() const
{
  if(!checked){throw wzFailure();}
}
void wzTensorGrid::add(wzIndex i, wzFloat f)
{
  int l; wzFloat t;
  if(checked){throw wzFailure();}
  if(f <= -wzInfty) return;
  if(f >=  wzInfty) return;
  wzRangeLoop(n[i],l){
    t=x[i][l];
    if(t < f-epsilon) continue;
    if(t < f+epsilon) return;
    x[i][l] = f; f = t;
  }
  l = n[i].create(); x[i][l] = f;
}

void wzTensorGrid::getBox(wzFloat *min, wzFloat *max) const
{
  for(int i=0;i<wzPointDim;i++){
    min[i] = x[i][n[i].first()];
    max[i] = x[i][n[i].last()];
  }
}

void wzTensorGrid::endInitialization(wzFloat *min, wzFloat *max)
{
  int i;
  if(checked){throw wzFailure();}
  for(i=wzPointDim-1;i>=0;i--){
    if(n[i].length()>1){
      dim = i+1; break;
    }
    if(n[i].length()==0) add(i,min[i]);
    if(n[i].length()==0) add(i,max[i]);
    if(n[i].length()==0) add(i,0);
  }
  for(i=0;i<dim;i++){
    if(n[i].length()==0){add(i,min[i]); add(i,max[i]);}
    if(n[i].length()==0){add(i,0); add(i,1);}
    if(n[i].length()==1){add(i,min[i]);}
    if(n[i].length()==1){add(i,max[i]);}
    if(n[i].length()==1){
      if(min[i]>0) 	add(i,x[i][0]+1);
      else		add(i,0);
    }
    if(n[i].length()==1){add(i,1);}
  }
  checked = 1;
}