dv_posf.cc

Go to the documentation of this file.

 
#include "dv_posf.h"
#include "domain.h"
#include <iostream>
#include <cstdlib>
#include <numeric> //accumulate
 
 
dv_posf::dv_posf(domain    *line,
                 const      string s,
                 const bool Lt,
                 const bool Lo) {
 
    domn          = line;
    var_name      = s;
    L_transported = Lt;
    L_output      = Lo;
    d             = vector<double>(domn->ngrdf, 0.0);
 
    double dx = domn->pram->domainLength / domn->ngrd;
    d.at(0) = domn->pram->xDomainCenter - 0.5*domn->pram->domainLength;
    for(int i=1; i<domn->ngrdf; i++)
        d.at(i) = d.at(i-1) + dx;
 
}
 
 
void dv_posf::splitCell(const int isplt,
                        const int nsplt,
                        const vector<double> &cellFaces) {
 
    d.insert( d.begin() + isplt+1, nsplt, 0.0 );
    for(int i=isplt+1, j=1; i<=isplt+nsplt; i++, j++)
        d.at(i) = cellFaces.at(j);
 
}
 
void dv_posf::merge2cells(const int    imrg,
                          const double m1,
                          const double m2,
                          const bool   LconstVolume) {
 
    d.erase(d.begin() + imrg+1);
    
    if(LconstVolume || domn->pram->bcType=="WALL") {   //todo: generalize this (works for constant density flows only, and not spatial (due to velocity)).
        domn->pos->merge2cells(imrg, m1, m2, LconstVolume);
        return;
    }
    
    double invC = 1.0/domn->pram->cCoord;
    double C    = domn->pram->cCoord;
 
    vector<double> dxc;
    domn->mesher->setGridDxc(domn, dxc, domn->pram->cCoord);
    dxc.at(imrg) = (m1+m2)/domn->rho->d.at(imrg);
    if(domn->pram->Lspatial)
        dxc.at(imrg) /= domn->uvel->d.at(imrg);
 
 
    domn->mesher->setGridFromDxc(dxc);       // does pos also --> pos is done for each merge and also at the end in meshManager::merge2cells
 
 
    //if(domn->pram->bcType=="OUTFLOW") {
    //
    //    double V2tot = accumulate(dxc.begin(), dxc.end(), 0.0);
    //    double dmb;
    //    d[0] = -pow(0.5*V2tot, invC);
    //    for(int ie=1, iw=0, i=0; ie<d.size(); ie++, iw++, i++) {
    //        if(d[iw] <= 0.0) {
    //            dmb = pow(abs(d[iw]),C) - dxc[i];
    //            if(dmb >=0) d[ie] = -pow( dmb, invC);
    //            else        d[ie] =  pow(-dmb, invC);
    //        }
    //        else {
    //            d[ie] = pow(pow(d[iw],C) + dxc[i], invC);
    //        }
    //    }
    //}
 
 
    //else if(domn->pram->bcType == "Wall_OUT") {
    //    d[0] = 0.0;
    //    for(int ie=1, iw=0, i=0; ie<domn->ngrdf; ie++, iw++, i++)
    //        d[ie] = pow(pow(d[iw],C) + dxc[i], invC);
    //}
 
 
    //else {
    //    cout << endl << "ERROR: dv_posf::merge2cells: not setup for given bcType" << endl;
    //    exit(0);
    //}
 
}
 
 
 
void dv_posf::setDvFromRegion(const int i1, const int i2){
 
    // note, we are owned by the eddyline, so domn is eddl, so to get domn data, use domn->domn
    const vector<double> &domn_data = domn->domn->varMap.find(var_name)->second->d;
 
    if(i2 >= i1)
        d.assign(domn_data.begin()+i1, domn_data.begin()+i2+2  );
    else {           // wrap around (periodic assignment)
        d.assign(domn_data.begin()+i1, domn_data.end()-1);
        double idmb  = d.size();
        d.insert(d.end(), domn_data.begin(), domn_data.begin()+i2+2 );
        for(int i=idmb; i<d.size(); i++)
            d.at(i)+=domn->domn->Ldomain();
    }
 
}
 
 
void dv_posf::resize() {
    d.resize(domn->ngrdf);
}