VM2D::VelocityBiotSavart Class Reference

Класс, определяющий способ вычисления скоростей More...

#include <Velocity2DBiotSavart.h>

Inheritance diagram for VM2D::VelocityBiotSavart:

Collaboration diagram for VM2D::VelocityBiotSavart:

Public Member Functions
	VelocityBiotSavart (const World2D &W_)
	Конструктор
virtual	~VelocityBiotSavart ()
	Деструктор
virtual float	CalcConvVeloToSetOfPointsFromWake (const WakeDataBase &pointsDb, std::vector< Point2D > &velo, std::vector< double > &domainRadius, bool calcVelo, bool calcRadius) override
	Вычисление конвективных скоростей и радиусов вихревых доменов в заданном наборе точек от следа
void	CalcConvVelo ()
	Вычисление конвективных скоростей вихрей и виртуальных вихрей в вихревом следе, а также в точках wakeVP.
virtual void	CalcConvVPVeloToSetOfPointsFromWake (const WakeDataBase &pointsDb, std::vector< Point2D > &velo, std::vector< double > &domainRadius, bool calcVelo, bool calcRadius)
virtual void	CalcVeloToWakeVP ()
	Вычисление скоростей в точках wakeVP.
void	CalcDiffVeloI1I2ToSetOfPointsFromWake (const WakeDataBase &pointsDb, const std::vector< double > &domainRadius, const WakeDataBase &vorticesDb, std::vector< double > &I1, std::vector< Point2D > &I2)
	Вычисление числителей и знаменателей диффузионных скоростей в заданном наборе точек
void	CalcDiffVeloI1I2ToSetOfPointsFromSheets (const WakeDataBase &pointsDb, const std::vector< double > &domainRadius, const Boundary &bnd, std::vector< double > &I1, std::vector< Point2D > &I2)
void	CalcDiffVeloI1I2ToWakeFromSheets (const WakeDataBase &pointsDb, const std::vector< double > &domainRadius, const Boundary &bnd, std::vector< double > &I1, std::vector< Point2D > &I2)
void	CalcDiffVeloI1I2ToWakeFromWake (const WakeDataBase &pointsDb, const std::vector< double > &domainRadius, const WakeDataBase &vorticesDb, std::vector< double > &I1, std::vector< Point2D > &I2)
void	CalcDiffVeloI1I2 ()
	Вычисление диффузионных скоростей вихрей и виртуальных вихрей в вихревом следе
void	CalcDiffVeloI0I3 ()
void	LimitDiffVelo (std::vector< Point2D > &diffVel)
	Контроль больших значений диффузионных скоростей
void	CalcDiffVelo ()
	Вычисление диффузионных скоростей
void	ResizeAndZero ()
	Очистка старых массивов под хранение скоростей, выделение новой памяти и обнуление
void	SaveVisStress ()
	Сохранение вязких напряжений
void	GetWakeInfluenceToRhs (const Airfoil &afl, std::vector< double > &wakeRhs) const
	Генерация вектора влияния вихревого следа на профиль
void	CPUGetFASTWakeInfluenceToRhs (std::vector< double > &wakeRhs, std::vector< double > &wakeRhsLin) const
void	FillRhs (Eigen::VectorXd &rhsReord) const

Public Attributes
VortexesParams	wakeVortexesParams
	Струтура, определяющая параметры вихрей в следе
std::vector< VortexesParams >	virtualVortexesParams
	Вектор струтур, определяющий параметры виртуальных вихрей для профилей
OptimizedVelocity	optimizedVelocity

Protected Attributes
const World2D &	W
	Константная ссылка на решаемую задачу

Detailed Description

Класс, определяющий способ вычисления скоростей

Способ вычисления скоростей

• напрямую по закону Био — Савара

Author

Марчевский Илья Константинович

Сокол Ксения Сергеевна

Рятина Евгения Павловна

Колганова Александра Олеговна

Version

1.14

Date

6 марта 2026 г.

Definition at line 64 of file Velocity2DBiotSavart.h.

Constructor & Destructor Documentation

VelocityBiotSavart()

VelocityBiotSavart::VelocityBiotSavart

(

const World2D &

W_

)

Конструктор

Parameters

[in]

W_

константная ссылка на решаемую задачу

Definition at line 54 of file Velocity2DBiotSavart.cpp.

                                                        :
    Velocity(W_)
{
};

~VelocityBiotSavart()

VelocityBiotSavart::~VelocityBiotSavart ( )

virtual

Деструктор

Definition at line 60 of file Velocity2DBiotSavart.cpp.

{   
};

Member Function Documentation

CalcConvVelo()

void Velocity::CalcConvVelo ( )

inherited

Вычисление конвективных скоростей вихрей и виртуальных вихрей в вихревом следе, а также в точках wakeVP.

Warning

скорости приплюсовываются к уже имеющимся

FMM ///

Definition at line 952 of file Velocity2D.cpp.

{
    W.getTimers().start("ConvVel");
    
    //Влияние следа на след
#if (defined(__CUDACC__) || defined(USE_CUDA)) && (defined(CU_CONV_TOWAKE))     
    //Вызов виртуальной функции
 
    W.timerConvVelo.reset();
    W.timerConvVelo.start();
 
    VMlib::vmTimer timerAgpu, timerBgpu;
 
    std::unique_ptr<BHcu::CudaTreeInfo>& cntrTree = W.getNonConstCuda().cntrTreeWake;
    auto& treeWake = *W.getNonConstCuda().inflTreeWake;
    if (W.getPassport().numericalSchemes.velocityComputation.second == 1)
    {
        cntrTree->MemoryAllocate((int)W.getCuda().n_CUDA_wake);
        cntrTree->Update((int)W.getWake().vtx.size(), W.getWake().devVtxPtr);
        cntrTree->Build();
 
        //Перестроение дерева для вихрей //todo //временно для корректного расчета epsast
        treeWake.MemoryAllocate((int)W.getCuda().n_CUDA_wake);
        treeWake.Update((int)W.getWake().vtx.size(), W.getWake().devVtxPtr);
        timerAgpu.start();      
        treeWake.Build();
        timerAgpu.stop();
 
        timerBgpu.start();
        treeWake.UpwardTraversal(W.getPassport().numericalSchemes.nbodyMultipoleOrder);
        timerBgpu.stop();
    }
        
    double timerCgpuduration = GPUCalcConvVeloToSetOfPointsFromWake(W.getNonConstCuda().cntrTreeWake, W.getWake(), wakeVortexesParams.convVelo, wakeVortexesParams.epsastWake, true, true);
    
 
    W.timerConvVelo.stop();
 
    /*
    std::vector<fmm::particle2d> sources(W.getWake().vtx.size());
    #pragma omp parallel for
    for (int i = 0; i < (int)sources.size(); ++i)
    {
        const auto& w = W.getWake().vtx[i];
    
        sources[i].center = fmm::point2d(w.r()[0], w.r()[1]);
        sources[i].q = w.g();
        sources[i].morton_code = 0;
    }
    
    auto targets = sources;
    for (int i = 0; i < (int)sources.size(); ++i)
        targets[i].q = 0;
 
    double t = omp_get_wtime();
    //fmm::ComputeExact<fmm::InteractionType::Potential2d>(sources, targets); // or fmm::gpu::ComputeExact... for gpu version
    fmm::ComputeExact<fmm::InteractionType::Force2d>(sources, targets); // or fmm::gpu::ComputeExact... for gpu version
    std::cout << "Compute exact time = " << omp_get_wtime() - t << std::endl;
    
    t = omp_get_wtime();    
 
    fmm::FastMultipole fmm(sources, targets, 1.e-8, fmm::FMM_AUTO, fmm::FMM_AUTO); // or fmm::gpu::FastMultipole... for gpu version
    fmm::ReadError<fmm::InteractionType::Force2d>(fmm.forces);
    //fmm::ReadError<fmm::InteractionType::Potential2d>(fmm.potentials);    
    std::cout << "Total time = " << omp_get_wtime() - t << std::endl;
    */
    
 
 
    /*
    std::ofstream treeTimeFile;
    if (W.getCurrentStep() == 0)
    {
        treeTimeFile.open(W.getPassport().dir + "/dbg/treeTime.csv");
        treeTimeFile << "step,time,N,lev,averPoints,tBUI,tUPW,tDNW,tBUIgpu,tUPWgpu,tDNWgpu\n";
    }
    else
        treeTimeFile.open(W.getPassport().dir + "/dbg/treeTime.csv", std::ios::app);
 
 
    int estLevel = std::max(4, (int)(log2(W.getWake().vtx.size())));
 
    for (int controlLevel = std::max(estLevel - 2, 1); controlLevel <= estLevel + 2; ++controlLevel)
    {
        VMlib::vmTimer timerA, timerB;
        if (W.getWake().vtx.size() > 0)
        {
            if (W.getPassport().numericalSchemes.velocityComputation.second == 1)
            {
                W.getInflTreeWake().Update(W.getWake().vtx, controlLevel);
                timerA.start();
                W.getInflTreeWake().Build();
                timerA.stop();
                timerB.start();
                W.getInflTreeWake().UpwardTraversal(W.getPassport().numericalSchemes.nbodyMultipoleOrder);
                timerB.stop();
 
                W.getCntrTreeWake().Update(W.getWake().vtx, controlLevel);
                W.getCntrTreeWake().Build();
                W.getCntrTreeWake().UpwardTraversal(W.getPassport().numericalSchemes.nbodyMultipoleOrder);
            }
 
            std::vector<Point2D> tempVelo(wakeVortexesParams.convVelo.size());
            std::vector<double> tempRadius(wakeVortexesParams.epsastWake.size());
 
            //Вызов виртуальной функции             
            double timerCduration = CalcConvVeloToSetOfPointsFromWake(W.getWake(), tempVelo, tempRadius, true, false);
 
            int nObjInContr = 0;
            for (int c = 0; c < W.getCntrTreeWake().indexControlCells.size(); ++c)
            {
                const auto& idx = W.getCntrTreeWake().indexControlCells[c];
                if (idx > W.getCntrTreeWake().object.size())
                    nObjInContr += 1;
                else
                {
                    nObjInContr += W.getCntrTreeWake().range[idx].second - W.getCntrTreeWake().range[idx].first + 1;
                }
            }
            
            treeTimeFile << W.getCurrentStep() << ',' << W.getCurrentTime() << ',' << W.getInflTreeWake().object.size() << ',' \
                << controlLevel << ',' << (double)nObjInContr / W.getCntrTreeWake().indexControlCells.size() << ',' \
                << timerA.duration() << ',' << timerB.duration() << ',' << timerCduration << ',' \
                << timerAgpu.duration() << ',' << timerBgpu.duration() << ',' << timerCgpuduration << '\n';
        }
    }
 
    treeTimeFile.close();
    */
#else
    float timeUpd1=0, timeBld1=0, timeUpw1=0, timeUpd2=0, timeBld2=0, timeUpw2=0, timeDnw=0;
    W.timerConvVelo.reset();
    W.timerConvVelo.start();
 
    VMlib::vmTimer fullTimer;
    fullTimer.reset();
    fullTimer.start();
 
    if (W.getWake().vtx.size() > 0)
    {
        if (W.getPassport().numericalSchemes.velocityComputation.second == 1)
        {
            int estOptLevel = std::max(4, (int)(log2(W.getWake().vtx.size())) - 2);
 
            timeUpd1 = W.getInflTreeWake().Update(W.getWake().vtx, estOptLevel);
            timeBld1 = W.getInflTreeWake().Build();
            timeUpw1 = W.getInflTreeWake().UpwardTraversal(W.getPassport().numericalSchemes.nbodyMultipoleOrder);
 
            //timeUpd2 = W.getCntrTreeWake().Update(W.getWake().vtx, estOptLevel);
            //timeBld2 = W.getCntrTreeWake().Build();
            //timeUpw2 = W.getCntrTreeWake().UpwardTraversal(W.getPassport().numericalSchemes.nbodyMultipoleOrder);
            
            //std::cout << "timeUpd = " << timeUpd1 << " " << timeUpd2 << '\n';
            //std::cout << "timeBld = " << timeBld1 << " " << timeBld2 << '\n';
            //std::cout << "timeUpw = " << timeUpw1 << " " << timeUpw2 << '\n';
        
        }
        
        //Вызов виртуальной функции 
        
        timeDnw = CalcConvVeloToSetOfPointsFromWake(W.getWake(), wakeVortexesParams.convVelo, wakeVortexesParams.epsastWake, true, true);
        
        //std::string vname = VMlib::fileNameStep("VeloBH", 0, wakeVortexesParams.convVelo.size(), "txt");
        //std::ofstream velofile;
        //velofile.open(W.getPassport().dir + vname);
 
        //for (int i = 0; i < wakeVortexesParams.convVelo.size(); ++i)
        //  velofile << wakeVortexesParams.convVelo[i][0] << " " << wakeVortexesParams.convVelo[i][1] << "\n";
        //velofile.close();
        //exit(-124);
 
        //std::cout << "timeDnw = " << timeDnw << '\n';
    }
    fullTimer.stop();
    
    W.timerConvVelo.stop();
    //std::cout << "theta = " << W.getPassport().numericalSchemes.nbodyTheta << "\n";
    
    /*
    VMlib::vmTimer timerExact;
 
    if (W.getWake().vtx.size() > 0)
    {
        std::vector<fmm::particle2d> sources(W.getWake().vtx.size());
        #pragma omp parallel for
        for (int i = 0; i < (int)sources.size(); ++i)
        {
            const auto& w = W.getWake().vtx[i];
 
            sources[i].center = fmm::point2d(w.r()[0], w.r()[1]);
            sources[i].q = w.g();
            sources[i].morton_code = 0;
        }
 
        auto targets = sources;
        for (int i = 0; i < (int)sources.size(); ++i)
            targets[i].q = 0;
 
        //double t = omp_get_wtime();
        //fmm::ComputeExact<fmm::InteractionType::Potential2d>(sources, targets); // or fmm::gpu::ComputeExact... for gpu version
        timerExact.start();
        fmm::ComputeExact<fmm::InteractionType::Force2d>(sources, targets); // or fmm::gpu::ComputeExact... for gpu version
        timerExact.stop();
        //std::cout << "Compute exact time = " << omp_get_wtime() - t << std::endl;
 
        //t = omp_get_wtime();
 
        fmm::FastMultipole fmm(sources, targets, 1.e-8, fmm::FMM_AUTO, fmm::FMM_AUTO); // or fmm::gpu::FastMultipole... for gpu version
        fmm::ReadError<fmm::InteractionType::Force2d>(fmm.forces);
        //fmm::ReadError<fmm::InteractionType::Potential2d>(fmm.potentials);
        //std::cout << "Total time = " << omp_get_wtime() - t << std::endl;
 
        std::cout << "Hybrid time ConvVelo Wake = " << fullTimer.duration() <<
            //W.timerConvVelo.duration() << 
            " ( " << timeUpd1 + timeUpd2 + timeBld1 + timeBld2 + timeUpw1 + timeUpw2 + timeDnw << " )" << '\n';
 
        fmm::ReadErrorHyb(*((std::vector<std::array<double, 2>>*) &wakeVortexesParams.convVelo));
        std::cout << "Exact time = " << timerExact.duration() << '\n';
    }
    //*/
 
 
 
#endif
 
    std::vector<Point2D> nullVector(0);
    //вычисление конвективных скоростей по закону Био-Савара от виртуальных вихрей (точнее, от слоев)
    for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
    {
        //Влияние на след от панелей
#if (defined(__CUDACC__) || defined(USE_CUDA)) && (defined(CU_CONVVIRT))    
 
        switch (W.getPassport().numericalSchemes.velocityComputation.second)
        {
        case 0:
            W.getBoundary(bou).GPUCalcConvVelocityToSetOfPointsFromSheets(W.getWake(), wakeVortexesParams.convVelo);
            break;
        case 1:
            if (bou == 0)
                GPUCalcConvVelocityToSetOfPointsFromSheets(W.getNonConstCuda().cntrTreeWake, W.getWake(), wakeVortexesParams.convVelo);
            break;
        }
        //      
 
#else
        W.getBoundary(bou).CalcConvVelocityToSetOfPointsFromSheets(W.getWake(), wakeVortexesParams.convVelo);
#endif
    }
 
    //Скорости только что рожденных вихрей
    {
        size_t nVirtVortices = 0;
        for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
            nVirtVortices += W.getBoundary(bou).virtualWake.vtx.size();
 
        size_t counter = wakeVortexesParams.convVelo.size() - nVirtVortices;
 
        for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
            for (size_t v = 0; v < W.getBoundary(bou).virtualWake.vtx.size(); ++v)
            {
                wakeVortexesParams.convVelo[counter] = W.getBoundary(bou).afl.getV(W.getBoundary(bou).virtualWake.aflPan[v].second) \
                    + W.getBoundary(bou).virtualWake.vecHalfGamma[v] \
                    - W.getV0();
                ++counter;
            }
    }
 
    //Заполнение структур данных виртуальных вихрей
    {
        size_t nVirtVortices = 0;
        for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
            nVirtVortices += W.getBoundary(bou).virtualWake.vtx.size();
 
        size_t counter = wakeVortexesParams.convVelo.size() - nVirtVortices;
 
        for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
            for (size_t v = 0; v < W.getBoundary(bou).virtualWake.vtx.size(); ++v)
            {
                virtualVortexesParams[bou].convVelo[v] = wakeVortexesParams.convVelo[counter];
                virtualVortexesParams[bou].epsastWake[v] = wakeVortexesParams.epsastWake[counter];
                ++counter;
            }
    }
 
    W.getTimers().stop("ConvVel");
    W.getTimers().start("VelPres");
 
    //вычисление скоростей в заданных точках только на соответствующем шаге по времени
    if ((W.getPassport().timeDiscretizationProperties.saveVPstep > 0) && (!(W.getCurrentStep() % W.getPassport().timeDiscretizationProperties.saveVPstep)))
#ifdef OPTIMIZER    
        if ((W.getCurrentStep() >= OPTIMIZER_START_STEP && W.getCurrentStep() <= OPTIMIZER_STOP_STEP) || W.getCurrentStep() == 0)
            CalcVeloToWakeVP();
#else
 
    //t1 = -omp_get_wtime();
    CalcVeloToWakeVP();
    //t1 += omp_get_wtime();
    //std::cout << "time ConvVelo VP = " << t1 << std::endl;
#endif
 
    /*if ((W.getPassport().timeDiscretizationProperties.saveVPstep > 0) && (!(W.getCurrentStep() % W.getPassport().timeDiscretizationProperties.saveVPstep)))
    {
        int currentStep = W.getCurrentStep();
        if (currentStep >= 1100 && currentStep <= 1300)
        {
            CalcVeloToWakeVP();
            optimizedVelocity.AddVelocities(W.getNonConstMeasureVP().getNonConstVelocity());
        }
 
        if (currentStep == 1200)
        {
            optimizedVelocity.PerformAveragingVelo(W.getPassport().dir);
        }
    }*/
 
    W.getTimers().stop("VelPres");
}//CalcConvVelo()

Here is the call graph for this function:

CalcConvVeloToSetOfPointsFromWake()

float VelocityBiotSavart::CalcConvVeloToSetOfPointsFromWake ( const WakeDataBase &  pointsDb, std::vector< Point2D > &  velo, std::vector< double > &  domainRadius, bool  calcVelo, bool  calcRadius  )

overridevirtual

Вычисление конвективных скоростей и радиусов вихревых доменов в заданном наборе точек от следа

Parameters

[in]	pointsDb	константная ссылка на базу данных пелены из вихрей, в которых надо сосчитать конвективные скорости
[out]	velo	ссылка на вектор скоростей в требуемых точках
[out]	domainRadius	ссылка на вектор радиусов вихревых доменов
[in]	calcVelo	признак вычисления скоростей в точках
[in]	calcRadius	признак вычисления радиусов доменов

TESTONLYVELO

Implements VM2D::Velocity.

Definition at line 92 of file Velocity2DBiotSavart.cpp.

{
    std::vector<Point2D> selfVelo(pointsDb.vtx.size());
    domainRadius.resize(pointsDb.vtx.size());
 
    double cft = IDPI;
 
#pragma warning (push)
#pragma warning (disable: 4101)
    //Локальные переменные для цикла
    Point2D velI;
    Point2D tempVel;
    double dst2eps, dst2;
#pragma warning (pop)
 
    VMlib::vmTimer timer;
    timer.start();
 
    if (calcVelo)
    {
#pragma omp parallel for default(none) shared(selfVelo, cft, calcVelo, calcRadius, pointsDb, domainRadius) private(tempVel, velI, dst2, dst2eps) schedule(dynamic, DYN_SCHEDULE)
        for (int i = 0; i < pointsDb.vtx.size(); ++i)
        {
            double ee2[3] = { 10000.0, 10000.0, 10000.0 };
 
            velI.toZero();
 
            const Point2D& posI = pointsDb.vtx[i].r();
 
            for (size_t j = 0; j < W.getWake().vtx.size(); ++j)
            {
                const Point2D& posJ = W.getWake().vtx[j].r();
 
                dst2 = (posI - posJ).length2();
 
                //Модифицируем массив квадратов расстояний до ближайших вихрей из wake
#ifndef TESTONLYVELO
                if (calcRadius)
                    VMlib::ModifyE2(ee2, dst2);
#endif 
 
                const double& gamJ = W.getWake().vtx[j].g();
                const double sigmaJ2 = sqr(W.getWake().vtx[j].sigma());
 
                tempVel.toZero();
                dst2eps = VMlib::boundDenom(dst2, sigmaJ2); //Сглаживать надо!!!                
 
                tempVel = { -posI[1] + posJ[1], posI[0] - posJ[0] };
                tempVel *= (gamJ / dst2eps);
                velI += tempVel;
            }
 
            for (size_t j = 0; j < W.getSource().vtx.size(); ++j)
            {
                const Point2D& posJ = W.getSource().vtx[j].r();
                const double& gamJ = W.getPassport().physicalProperties.accelCft(W.getCurrentTime()) * W.getSource().vtx[j].g();
                const double sigmaJ2 = sqr(W.getSource().vtx[j].sigma());
 
                tempVel.toZero();
 
                dst2 = dist2(posI, posJ);
                dst2eps = VMlib::boundDenom(dst2, sigmaJ2); //Сглаживать надо!!!
 
                tempVel = { posI[0] - posJ[0], posI[1] - posJ[1] };
                tempVel *= (gamJ / dst2eps);
                velI += tempVel;
            }
#ifndef TESTONLYVELO
            if (calcRadius)
            {
                for (size_t s = 0; s < W.getNumberOfBoundary(); ++s)
                {
                    const auto& bou = W.getBoundary(s);
                    //Модифицируем массив квадратов расстояний до ближайших вихрей из virtualWake                   
                    for (size_t j = 0; j < bou.virtualWake.vtx.size(); ++j)
                        ModifyE2(ee2, dist2(posI, bou.virtualWake.vtx[j].r()));
                }
            }
#endif
 
            velI *= cft;
            selfVelo[i] = velI;
 
#ifndef TESTONLYVELO
            if (calcRadius)
                domainRadius[i] = 1.0 * sqrt((ee2[0] + ee2[1] + ee2[2]) / 3.0);
#endif
        }
    }
    else if (calcRadius)
    {
#pragma omp parallel for default(none) shared(selfVelo, cft, calcVelo, calcRadius, pointsDb, domainRadius) private(tempVel, velI, dst2) schedule(dynamic, DYN_SCHEDULE)
        for (int i = 0; i < pointsDb.vtx.size(); ++i)
        {
            double ee2[3] = { 10000.0, 10000.0, 10000.0 };
 
            velI.toZero();
 
            const Point2D& posI = pointsDb.vtx[i].r();
 
            for (size_t j = 0; j < W.getWake().vtx.size(); ++j)
            {
                const Point2D& posJ = W.getWake().vtx[j].r();
 
                dst2 = dist2(posI, posJ);
 
                //Модифицируем массив квадратов расстояний до ближайших вихрей из wake
                VMlib::ModifyE2(ee2, dst2);
            }
 
            for (size_t s = 0; s < W.getNumberOfBoundary(); ++s)
            {
                for (size_t j = 0; j < W.getBoundary(s).virtualWake.vtx.size(); ++j)
                {
                    const Point2D& posJ = W.getBoundary(s).virtualWake.vtx[j].r();
                    dst2 = dist2(posI, posJ);
 
                    //Модифицируем массив квадратов расстояний до ближайших вихрей из virtualWake
                    ModifyE2(ee2, dst2);
                }
            }
 
            domainRadius[i] = 1.0 * sqrt((ee2[0] + ee2[1] + ee2[2]) / 3.0);
        }
    } //else
 
 
    timer.stop();
 
    if (calcVelo)
        for (size_t i = 0; i < velo.size(); ++i)
            velo[i] += selfVelo[i];
 
#ifdef USE_CUDA 
    if (calcVelo)
        W.getCuda().CopyMemToDev<double, 2>(velo.size(), (double*)velo.data(), pointsDb.devVelPtr);
    if (calcRadius)
        W.getCuda().CopyMemToDev<double, 1>(domainRadius.size(), domainRadius.data(), pointsDb.devRadPtr);
#endif
    return (float)timer.duration();
}//CalcConvVeloToSetOfPointsFromWake(...)

Here is the call graph for this function:

CalcConvVPVeloToSetOfPointsFromWake()

virtual void VM2D::Velocity::CalcConvVPVeloToSetOfPointsFromWake ( const WakeDataBase &  pointsDb, std::vector< Point2D > &  velo, std::vector< double > &  domainRadius, bool  calcVelo, bool  calcRadius  )

inlinevirtualinherited

Reimplemented in VM2D::VelocityBarnesHut.

Definition at line 141 of file Velocity2D.h.

{};

CalcDiffVelo()

void Velocity::CalcDiffVelo ( )

inherited

Вычисление диффузионных скоростей

Вызывается в CalcVortexVelo()

Definition at line 155 of file Velocity2D.cpp.

{
    W.getTimers().start("DiffVel");
 
    double t12start, t12finish;
    double t03start, t03finish;
    double tOtherstart, tOtherfinish;
 
    if (W.getPassport().physicalProperties.nu > 0.0)
    {
        t12start = omp_get_wtime();
        CalcDiffVeloI1I2();
        t12finish = omp_get_wtime();
 
        t03start = omp_get_wtime();
        CalcDiffVeloI0I3();
        t03finish = omp_get_wtime();
  
        tOtherstart = omp_get_wtime();
        for (Point2D& diffV : wakeVortexesParams.diffVelo)
            diffV *= W.getPassport().physicalProperties.nu;
 
        //контроль застрелов диффузионной скорости      
        LimitDiffVelo(wakeVortexesParams.diffVelo);
 
        for (size_t afl = 0; afl < W.getNumberOfAirfoil(); ++afl)
            for (size_t i = 0; i < W.getAirfoil(afl).viscousStress.size(); ++i)
                W.getNonConstAirfoil(afl).viscousStress[i] *= W.getPassport().physicalProperties.nu;
 
        //SaveVisStress();
        
 
        //Заполнение структуры данных виртуальных вихрей
        size_t nVirtVortices = 0;
        for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
            nVirtVortices += W.getBoundary(bou).virtualWake.vtx.size();
 
        size_t counter = wakeVortexesParams.convVelo.size() - nVirtVortices;
 
        for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
            for (size_t v = 0; v < W.getBoundary(bou).virtualWake.vtx.size(); ++v)
            {           
                virtualVortexesParams[bou].diffVelo[v] = wakeVortexesParams.diffVelo[counter];
                virtualVortexesParams[bou].I0[v] = wakeVortexesParams.I0[counter];
                virtualVortexesParams[bou].I1[v] = wakeVortexesParams.I1[counter];
                virtualVortexesParams[bou].I2[v] = wakeVortexesParams.I2[counter];
                virtualVortexesParams[bou].I3[v] = wakeVortexesParams.I3[counter];              
                ++counter;
            }
 
        tOtherfinish = omp_get_wtime();     
    }
    W.getTimers().stop("DiffVel");
 
}// CalcDiffVelo()

Here is the call graph for this function:

CalcDiffVeloI0I3()

void Velocity::CalcDiffVeloI0I3 ( )

inherited

Definition at line 106 of file Velocity2D.cpp.

{
    for (size_t afl = 0; afl < W.getNumberOfAirfoil(); ++afl)
    {
#if (defined(__CUDACC__) || defined(USE_CUDA)) && (defined(CU_I0I3))
        W.getNonConstAirfoil(afl).GPUGetDiffVelocityI0I3ToSetOfPointsAndViscousStresses(wakeVortexesParams.epsastWake, wakeVortexesParams.I0, wakeVortexesParams.I3);
#else
        W.getNonConstAirfoil(afl).GetDiffVelocityI0I3ToSetOfPointsAndViscousStresses(W.getWake(), wakeVortexesParams.epsastWake, wakeVortexesParams.I0, wakeVortexesParams.I3);
#endif
    }   
 
 
    //влияние поверхности
    Point2D I3;
    double I0;
 
    double domrad = 0.0;
 
#pragma omp parallel for private(I0, I3, domrad)
    for (int vt = 0; vt < (int)wakeVortexesParams.diffVelo.size(); ++vt)
    {
        domrad = std::max(wakeVortexesParams.epsastWake[vt], W.getPassport().wakeDiscretizationProperties.getMinEpsAst());
 
        wakeVortexesParams.I0[vt] *= domrad;
        wakeVortexesParams.I0[vt] += DPI * sqr(domrad);
            
 
        I3 = wakeVortexesParams.I3[vt];
        I0 = wakeVortexesParams.I0[vt];
        
 
        if (fabs(I0) > 1.e-8)
            wakeVortexesParams.diffVelo[vt] += I3 * (1.0 / I0);     
    }
}//CalcDiffVeloI0I3()

Here is the call graph for this function:

Here is the caller graph for this function:

CalcDiffVeloI1I2()

void Velocity::CalcDiffVeloI1I2 ( )

inherited

Вычисление диффузионных скоростей вихрей и виртуальных вихрей в вихревом следе

Вызывает 4 раза функцию CalcDiffVeloToSetOfPoints

Warning

скорости приплюсовываются к уже имеющимся

Definition at line 62 of file Velocity2D.cpp.

{   
    // !!! пелена на пелену
#if (defined(__CUDACC__) || defined(USE_CUDA)) && (defined(CU_I1I2))    
    if (W.getPassport().numericalSchemes.velocityComputation.second == 0)
        GPUCalcDiffVeloI1I2ToSetOfPointsFromWake(W.getWake(), wakeVortexesParams.epsastWake, W.getWake(), wakeVortexesParams.I1, wakeVortexesParams.I2);        
    else
        GPUDiffVeloFAST(wakeVortexesParams.epsastWake, wakeVortexesParams.I1, wakeVortexesParams.I2);       
#else
    CalcDiffVeloI1I2ToWakeFromWake(W.getWake(), wakeVortexesParams.epsastWake, W.getWake(), wakeVortexesParams.I1, wakeVortexesParams.I2);  
#endif
 
     
    for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
    {
        //не нужно, т.к. сделано выше перед началом вычисления скоростей
        //W.getNonConstBoundary(bou).virtualWake.WakeSynchronize();
        
        
    //виртуальные на границе на след
#if (defined(__CUDACC__) || defined(USE_CUDA)) && (defined(CU_I1I2))                
        GPUCalcDiffVeloI1I2ToSetOfPointsFromSheets(W.getWake(), wakeVortexesParams.epsastWake, W.getBoundary(bou), wakeVortexesParams.I1, wakeVortexesParams.I2);
#else           
        CalcDiffVeloI1I2ToWakeFromSheets(W.getWake(), wakeVortexesParams.epsastWake, W.getBoundary(bou), wakeVortexesParams.I1, wakeVortexesParams.I2);
#endif          
    } //for bou
    
    Point2D I2;
    double I1;
    
#pragma omp parallel for private(I1, I2)
    for (int vt = 0; vt < (int)wakeVortexesParams.diffVelo.size(); ++vt)
    {
        I2 = wakeVortexesParams.I2[vt];
        I1 = wakeVortexesParams.I1[vt];
 
        if (fabs(I1) < 1.e-8)
            wakeVortexesParams.diffVelo[vt] = { 0.0, 0.0 };
        else
            wakeVortexesParams.diffVelo[vt] = I2 * (1.0 / (I1 * std::max(wakeVortexesParams.epsastWake[vt], W.getPassport().wakeDiscretizationProperties.getMinEpsAst())));
    }   
}//CalcDiffVeloI1I2()

Here is the call graph for this function:

Here is the caller graph for this function:

CalcDiffVeloI1I2ToSetOfPointsFromSheets()

void Velocity::CalcDiffVeloI1I2ToSetOfPointsFromSheets ( const WakeDataBase &  pointsDb, const std::vector< double > &  domainRadius, const Boundary &  bnd, std::vector< double > &  I1, std::vector< Point2D > &  I2  )

inherited

Todo:

Понять природу магической константы 8.0 и синхронизировать с GPU

Todo:

Сделать переменной и синхронизировать с GPU

Todo:

Учитываем пока только нулевой момент решения

Definition at line 271 of file Velocity2D.cpp.

{
    double tCPUSTART, tCPUEND;
 
    tCPUSTART = omp_get_wtime();
 
    std::vector<double> selfI1(pointsDb.vtx.size(), 0.0);
    std::vector<Point2D> selfI2(pointsDb.vtx.size(), { 0.0, 0.0 });
 
#pragma warning (push)
#pragma warning (disable: 4101)
    //Локальные переменные для цикла
    Point2D Rij;
    double rij, expr;
    double diffRadius;
    double left;
    double right;
    double posJx;
    double domRad;
#pragma warning (pop)
 
 
#pragma omp parallel for default(none) shared(selfI1, selfI2, domainRadius, bnd, pointsDb, std::cout) private(Rij, rij, expr, domRad, diffRadius, left, right, posJx)
    for (int i = 0; i < pointsDb.vtx.size(); ++i)
    {
        const Vortex2D& vtxI = pointsDb.vtx[i];
 
        domRad = std::max(domainRadius[i], W.getPassport().wakeDiscretizationProperties.getMinEpsAst());
 
        diffRadius = 8.0 * domRad;
 
        left = vtxI.r()[0] - diffRadius;
        right = vtxI.r()[0] + diffRadius;
 
        for (size_t j = 0; j < bnd.afl.getNumberOfPanels(); ++j)
        {
            const int nQuadPt = 3;
 
            const double ptG = bnd.sheets.freeVortexSheet(j, 0) * bnd.afl.len[j] / nQuadPt;
 
            for (int q = 0; q < nQuadPt; ++q)
            {
                const Point2D& ptJ = bnd.afl.getR(j) + bnd.afl.tau[j] * (q + 0.5) * bnd.afl.len[j] * (1.0 / nQuadPt);  // vorticesDb.vtx[j];
                posJx = ptJ[0];
 
                if ((left < posJx) && (posJx < right))
                {
                    Rij = vtxI.r() - ptJ;
                    rij = Rij.length();
                    if (rij < diffRadius && rij > 1.e-10)
                    {
                        expr = exp(-rij / domRad);
                        selfI2[i] += (ptG * expr / rij) * Rij;
                        selfI1[i] += ptG * expr;
                    }
                }//if (rij>1e-6)
            }
        }//for j
    } // for r
 
    for (size_t i = 0; i < I1.size(); ++i)
    {
        I1[i] += selfI1[i];
        I2[i] += selfI2[i];
    }
 
    tCPUEND = omp_get_wtime();
    //W.getInfo('t') << "DIFF_CPU: " << tCPUEND - tCPUSTART << std::endl;
}//CalcDiffVeloI1I2ToSetOfPointsFromSheets(...)

Here is the call graph for this function:

Here is the caller graph for this function:

CalcDiffVeloI1I2ToSetOfPointsFromWake()

void Velocity::CalcDiffVeloI1I2ToSetOfPointsFromWake ( const WakeDataBase &  pointsDb, const std::vector< double > &  domainRadius, const WakeDataBase &  vorticesDb, std::vector< double > &  I1, std::vector< Point2D > &  I2  )

inherited

Вычисление числителей и знаменателей диффузионных скоростей в заданном наборе точек

Parameters

[in]	pointsDb	константная ссылка на базу данных пелены из вихрей, в которых надо сосчитать диффузионные скорости
[in]	domainRadius	константная ссылка на вектор радиусов вихревых доменов
[in]	vorticesDb	константная ссылка на на базу данных пелены из вихрей,от которых надо сосчитать влияния на points
[out]	I1	ссылка на вектор величин I1 (знаменателей в диффузионных скоростях) в требуемых точках
[out]	I2	ссылка на вектор величин I2 (числителей в диффузионных скоростях) в требуемых точках

Todo:

Понять природу магической константы 8.0 и синхронизировать с GPU

Definition at line 213 of file Velocity2D.cpp.

{
    std::vector<double> selfI1(pointsDb.vtx.size(), 0.0);
    std::vector<Point2D> selfI2(pointsDb.vtx.size(), { 0.0, 0.0 }); 
 
#pragma warning (push)
#pragma warning (disable: 4101)
    //Локальные переменные для цикла
    Point2D Rij;
    double rij, expr;
    double diffRadius, domRad;
    double left;
    double right;
    double posJx;
#pragma warning (pop)
 
#pragma omp parallel for default(none) shared(selfI1, selfI2, domainRadius, vorticesDb, pointsDb) private(Rij, rij, expr, diffRadius, domRad, left, right, posJx)
    for (int i = 0; i < pointsDb.vtx.size(); ++i)
    {
        const Vortex2D& vtxI = pointsDb.vtx[i];
 
        domRad = std::max(domainRadius[i], W.getPassport().wakeDiscretizationProperties.getMinEpsAst());
 
        diffRadius = 8.0 * domRad;
 
        left = vtxI.r()[0] - diffRadius;
        right = vtxI.r()[0] + diffRadius;
 
        for (size_t j = 0; j < vorticesDb.vtx.size(); ++j)
        {
            const Vortex2D& vtxJ = vorticesDb.vtx[j];
            posJx = vtxJ.r()[0];
 
            if ((left < posJx) && (posJx < right))
            {
                Rij = vtxI.r() - vtxJ.r();
                rij = Rij.length();
                if (rij < diffRadius && rij > 1.e-10)
                {
                    expr = exp(-rij / domRad);
                    selfI2[i] += (vtxJ.g()* expr / rij) * Rij;
                    selfI1[i] += vtxJ.g()*expr;
                }
            }//if (rij>1e-6)
        }//for j
    } // for r
 
 
    for (size_t i = 0; i < I1.size(); ++i)
    {
        I1[i] += selfI1[i];
        I2[i] += selfI2[i];
    }
}//CalcDiffVeloI1I2ToSetOfPointsFromWake(...)

Here is the call graph for this function:

Here is the caller graph for this function:

CalcDiffVeloI1I2ToWakeFromSheets()

void Velocity::CalcDiffVeloI1I2ToWakeFromSheets ( const WakeDataBase &  pointsDb, const std::vector< double > &  domainRadius, const Boundary &  bnd, std::vector< double > &  I1, std::vector< Point2D > &  I2  )

inherited

Definition at line 1372 of file Velocity2D.cpp.

{
    CalcDiffVeloI1I2ToSetOfPointsFromSheets(pointsDb, domainRadius, bnd, I1, I2);
}

Here is the call graph for this function:

Here is the caller graph for this function:

CalcDiffVeloI1I2ToWakeFromWake()

void Velocity::CalcDiffVeloI1I2ToWakeFromWake ( const WakeDataBase &  pointsDb, const std::vector< double > &  domainRadius, const WakeDataBase &  vorticesDb, std::vector< double > &  I1, std::vector< Point2D > &  I2  )

inherited

Definition at line 1367 of file Velocity2D.cpp.

{
    CalcDiffVeloI1I2ToSetOfPointsFromWake(pointsDb, domainRadius, vorticesDb, I1, I2);
}//CalcDiffVeloI1I2ToWakeFromWake(...)

Here is the call graph for this function:

Here is the caller graph for this function:

CalcVeloToWakeVP()

void Velocity::CalcVeloToWakeVP ( )

virtualinherited

Вычисление скоростей в точках wakeVP.

Definition at line 1278 of file Velocity2D.cpp.

{
    std::vector<Point2D> velConvWake;
    std::vector<std::vector<Point2D>> velConvBou;
 
    int addWSize = (int)W.getMeasureVP().getWakeVP().vtx.size();
 
    velConvWake.resize(addWSize, { 0.0, 0.0 });
 
    velConvBou.resize(W.getNumberOfBoundary());
    for (size_t i = 0; i < W.getNumberOfBoundary(); ++i)
        velConvBou[i].resize(addWSize, { 0.0, 0.0 });
 
    std::vector<Point2D>& velocityRef = W.getNonConstMeasureVP().getNonConstVelocity();
    velocityRef.assign(addWSize, W.getV0());
 
#if (defined(USE_CUDA))
    W.getNonConstCuda().RefreshVP();
#endif  
 
#if (defined(__CUDACC__) || defined(USE_CUDA)) && (defined(CU_VP))
    std::unique_ptr<BHcu::CudaTreeInfo>& cntrTree = W.getNonConstCuda().cntrTreeVP;
    
    switch(W.getPassport().numericalSchemes.velocityComputation.second)
    {
    case 1:
        if (W.getCurrentStep() == 0 || W.isAnyMovableOrDeformable())
        {           
            cntrTree->MemoryAllocate((int)W.getCuda().n_CUDA_velVP);            
            cntrTree->Update((int)W.getMeasureVP().getWakeVP().vtx.size(), W.getMeasureVP().getWakeVP().devVtxPtr);
            cntrTree->Build();  
            GPUCalcConvVeloToSetOfPointsFromWake(cntrTree, W.getMeasureVP().getWakeVP(), velocityRef, W.getNonConstMeasureVP().getNonConstDomainRadius(), true, false);
        }
        break;
    case 0:
        GPUCalcConvVeloToSetOfPointsFromWake(cntrTree, W.getMeasureVP().getWakeVP(), velocityRef, W.getNonConstMeasureVP().getNonConstDomainRadius(), true, false);
        break;
    }
    
 
#else
    switch (W.getPassport().numericalSchemes.velocityComputation.second)
    {
    case 0:
        CalcConvVeloToSetOfPointsFromWake(W.getMeasureVP().getWakeVP(), velConvWake, W.getNonConstMeasureVP().getNonConstDomainRadius(), true, false);
        break;
    case 1:
        int order = W.getPassport().numericalSchemes.nbodyMultipoleOrder;
        double theta = W.getPassport().numericalSchemes.nbodyTheta;
 
        W.getCntrTreeVP().Update(W.getMeasureVP().getWakeVP().vtx);
        W.getCntrTreeVP().Build();
        W.getCntrTreeVP().UpwardTraversal(W.getPassport().numericalSchemes.nbodyMultipoleOrder);
        W.getInflTreeWake().DownwardTraversalVorticesToPoints(W.getCntrTreeVP(), velConvWake, W.getNonConstMeasureVP().getNonConstDomainRadius(), theta, order, false);
    }
#endif  
 
     
    for (size_t i = 0; i < W.getNumberOfBoundary(); ++i)
#if (defined(__CUDACC__) || defined(USE_CUDA)) && (defined(CU_VP))          
        switch (W.getPassport().numericalSchemes.velocityComputation.second)
        {
        case 0:
            W.getNonConstBoundary(i).GPUCalcConvVelocityToSetOfPointsFromSheets(W.getMeasureVP().getWakeVP(), velConvBou[i]);
 
            for (int s = 0; s < addWSize; ++s)
                velocityRef[s] += velConvWake[s];
 
            for (size_t bou = 0; bou < velConvBou.size(); ++bou)
                for (int j = 0; j < addWSize; ++j)
                    velocityRef[j] += velConvBou[bou][j];
            break;
        case 1:
            if (i == 0)
                GPUCalcConvVelocityToSetOfPointsFromSheets(cntrTree, W.getMeasureVP().getWakeVP(), velocityRef);
            break;
        }
#else
        W.getNonConstBoundary(i).CalcConvVelocityToSetOfPointsFromSheets(W.getMeasureVP().getWakeVP(), velConvBou[i]);
        for (int i = 0; i < addWSize; ++i)
            velocityRef[i] += velConvWake[i];
 
        for (size_t bou = 0; bou < velConvBou.size(); ++bou)
            for (int j = 0; j < addWSize; ++j)
                velocityRef[j] += velConvBou[bou][j];
#endif
}// CalcVeloToWakeVP()

Here is the call graph for this function:

Here is the caller graph for this function:

CPUGetFASTWakeInfluenceToRhs()

void Velocity::CPUGetFASTWakeInfluenceToRhs ( std::vector< double > &  wakeRhs, std::vector< double > &  wakeRhsLin  ) const

inherited

Definition at line 656 of file Velocity2D.cpp.

{
    const size_t nvt = W.getWake().vtx.size();
 
 
    size_t nTotPan = 0;
 
    for (size_t s = 0; s < W.getNumberOfAirfoil(); ++s)
        nTotPan += W.getAirfoil(s).getNumberOfPanels();
 
    wakeRhs.assign(nTotPan, 0.0);
    wakeRhsLin.assign(nTotPan, 0.0);
 
    if (nvt > 0)
    {
        auto& inflTree = W.getInflTreeWake();
 
        auto& cntrTree = W.getCntrTreePnl();
 
        inflTree.DownwardTraversalVorticesToPanels(
            cntrTree,
            wakeRhs,
            wakeRhsLin,
            W.getPassport().numericalSchemes.nbodyTheta,
            W.getPassport().numericalSchemes.nbodyMultipoleOrder);
    }
}//CPUGetFASTWakeInfluenceToRhs(...)

Here is the call graph for this function:

Here is the caller graph for this function:

FillRhs()

void Velocity::FillRhs ( Eigen::VectorXd &  rhsReord ) const

inherited

!!!!!!!!!!влияние присоединенных слоев от самого себя и от других профилей

Definition at line 811 of file Velocity2D.cpp.

{
    Eigen::VectorXd locRhs;
    std::vector<double> lastRhs(W.getNumberOfBoundary());
 
    size_t currentRow = 0;
    size_t currentSkosRow = 0;
 
    size_t nAllVars = 0;
    for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
        nAllVars += W.getBoundary(bou).GetUnknownsSize();
 
    //Временные массивы для хранения правой части по быстрому методу (размеры задаются внутри CPUGetFASTWakeInfluenceToRhs)
    std::vector<double> fastWakeRhs;
    std::vector<double> fastWakeRhsLin;
 
#ifndef USE_CUDA
    if (W.getPassport().numericalSchemes.velocityComputation.second == 1) //Быстрый метод на CPU
        CPUGetFASTWakeInfluenceToRhs(fastWakeRhs, fastWakeRhsLin);  
#endif
    size_t curGlobPnl = 0;
 
    for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
    {
        const Airfoil& afl = W.getAirfoil(bou);
        size_t np = afl.getNumberOfPanels();
 
        size_t nVars;
 
        nVars = W.getBoundary(bou).GetUnknownsSize();
        locRhs.resize(nVars);
 
        std::vector<double> wakeRhs;
 
        double tt1 = omp_get_wtime();
 
#if (defined(__CUDACC__) || defined(USE_CUDA)) && (defined(CU_RHS))
        W.timerRhs.reset();
        W.timerRhs.start();
        if (W.getPassport().numericalSchemes.velocityComputation.second == 1)       
            GPUFASTGetWakeInfluenceToRhs(afl, wakeRhs);
        else
            GPUGetWakeInfluenceToRhs(afl, wakeRhs);
        W.timerRhs.stop();
#else
        if (W.getPassport().numericalSchemes.velocityComputation.second == 1) //Собираем результат из быстрого метода
        {
            size_t shDim = W.getBoundary(afl.numberInPassport).sheetDim;
 
            wakeRhs.assign(np * shDim, 0.0);
 
            for (size_t i = 0; i < np; ++i)
            {
                wakeRhs[i] = fastWakeRhs[curGlobPnl + i];
 
                if (shDim != 1)
                    wakeRhs[np + i] = fastWakeRhsLin[curGlobPnl + i];
            }
            curGlobPnl += np;
        }
        else
            GetWakeInfluenceToRhs(afl, wakeRhs); //Прямой расчет
 
        //std::string vname = VMlib::fileNameStep("rhsLinBH", 0, wakeRhs.size(), "txt");
        //std::ofstream velofile;
        //velofile.open(W.getPassport().dir + vname);
 
        //for (int i = 0; i < wakeRhs.size(); ++i)
        //  velofile << wakeRhs[i] << "\n";
        //velofile.close();
        //exit(-124);
#endif  
        std::vector<double> vInfRhs;
        afl.GetInfluenceFromVInfToPanel(vInfRhs);
 
#pragma omp parallel for \
    default(none) \
    shared(locRhs, afl, bou, wakeRhs, vInfRhs, np) \
    schedule(dynamic, DYN_SCHEDULE)
        for (int i = 0; i < (int)afl.getNumberOfPanels(); ++i)
        {
            locRhs(i) = -vInfRhs[i] - wakeRhs[i] + 0.25 * ((afl.getV(i) + afl.getV(i + 1)) & afl.tau[i]); //0.25 * (afl.getV(i) + afl.getV(i + 1))*afl.tau[i] - прямолинейные
            if (W.getBoundary(bou).sheetDim > 1)
                locRhs(np + i) = -vInfRhs[np + i] - wakeRhs[np + i];
 
            //Если включен быстрый метод и есть подвижные тела, то пока что все равно вычисляем IQ
            if (W.isAnyMovableOrDeformable())
            {
                for (size_t q = 0; q < W.getNumberOfBoundary(); ++q)
                {
                    const auto& sht = W.getBoundary(q).sheets;
                    const auto& iq = W.getIQ(bou, q);
 
                    const Airfoil& aflOther = W.getAirfoil(q);
                    if (W.getMechanics(q).isMoves)
                    {
                        for (size_t j = 0; j < aflOther.getNumberOfPanels(); ++j)
                        {
                            if ((i != j) || (bou != q))
                            {
                                locRhs(i) += -iq.first(i, j) * sht.attachedVortexSheet(j, 0);
                                locRhs(i) += -iq.second(i, j) * sht.attachedSourceSheet(j, 0);
                            }//if (i != j)
 
                        }//for j
                    }
                }//for q
            }
        }//for i
 
        lastRhs[bou] = 0.0;
 
#pragma omp for
        for (int q = 0; q < (int)afl.gammaThrough.size(); ++q)
            lastRhs[bou] += afl.gammaThrough[q];
 
 
        const double currT = W.getCurrentTime();
        const double dt = W.getPassport().timeDiscretizationProperties.dt;
 
        lastRhs[bou] += (W.getMechanics(bou).circulation - W.getMechanics(bou).circulationOld) * (W.getAirfoil(bou).inverse ? 1.0 : -1.0);
 
 
        //размазываем правую часть      
        for (size_t i = 0; i < nVars; ++i)
            rhsReord(i + currentSkosRow) = locRhs(i);
 
        rhsReord(nAllVars + bou) = lastRhs[bou];
 
        currentRow += nVars + 1;
        currentSkosRow += nVars;
    }// for bou
}

Here is the call graph for this function:

GetWakeInfluenceToRhs()

void Velocity::GetWakeInfluenceToRhs ( const Airfoil &  afl, std::vector< double > &  wakeRhs  ) const

inherited

Генерация вектора влияния вихревого следа на профиль

Генерирует вектор влияния вихревого следа на профиль, используемый затем для расчета вектора правой части.

Parameters

[in]	afl	константная ссылка на профиль, правая часть для которого вычисляется
[out]	wakeRhs	ссылка на вектор влияния вихревого следа на ОДИН профиль

Definition at line 619 of file Velocity2D.cpp.

{
    size_t np = afl.getNumberOfPanels();
    size_t shDim = W.getBoundary(afl.numberInPassport).sheetDim;
 
    wakeRhs.resize(W.getBoundary(afl.numberInPassport).GetUnknownsSize());
 
    //локальные переменные для цикла    
    std::vector<double> velI(shDim, 0.0);
 
#pragma omp parallel for default(none) shared(shDim, afl, np, wakeRhs, IDPI) private(velI)
    for (int i = 0; i < np; ++i)
    {
        velI.assign(shDim, 0.0);
 
        if (W.getWake().vtx.size() > 0)
        {
            //Учет влияния следа
            afl.GetInfluenceFromVorticesToPanel(i, W.getWake().vtx.data(), W.getWake().vtx.size(), velI);
        }
 
        if (W.getSource().vtx.size() > 0)
        {
            //Учет влияния источников
            afl.GetInfluenceFromSourcesToPanel(i, W.getSource().vtx.data(), W.getSource().vtx.size(), velI);
        }
 
        for (size_t j = 0; j < shDim; ++j)
            velI[j] *= IDPI / afl.len[i];
 
        wakeRhs[i] = velI[0];
 
        if (shDim != 1)
            wakeRhs[np + i] = velI[1];
    }//for i
}//GetWakeInfluenceToRhs(...)

Here is the call graph for this function:

Here is the caller graph for this function:

LimitDiffVelo()

void Velocity::LimitDiffVelo ( std::vector< Point2D > &  diffVel )

inherited

Контроль больших значений диффузионных скоростей

Parameters

[in,out]

diffVel

ссылка на вектор диффузионных скоростей

Definition at line 143 of file Velocity2D.cpp.

{
    for (size_t i = 0; i < diffVel.size(); ++i)
    {
        Point2D& diffV = diffVel[i];        
 
        if (diffV.length() > 1.5 * W.getPassport().physicalProperties.vRef)
            diffV.normalize(1.5 * W.getPassport().physicalProperties.vRef);
    }
}

Here is the call graph for this function:

Here is the caller graph for this function:

ResizeAndZero()

void Velocity::ResizeAndZero ( )

inherited

Очистка старых массивов под хранение скоростей, выделение новой памяти и обнуление

Вызывается в CalcVortexVelo() на каждом шаге расчета перед непосредственно расчетом скоростей

Definition at line 523 of file Velocity2D.cpp.

{
    W.getTimers().start("ConvVel"); 
    wakeVortexesParams.convVelo.clear();
    wakeVortexesParams.convVelo.resize(W.getWake().vtx.size(), { 0.0, 0.0 });
    W.getTimers().stop("ConvVel");
 
    W.getTimers().start("DiffVel");
    wakeVortexesParams.I0.clear();
    wakeVortexesParams.I0.resize(W.getWake().vtx.size(), 0.0);
 
    wakeVortexesParams.I1.clear();
    wakeVortexesParams.I1.resize(W.getWake().vtx.size(), 0.0);
 
    wakeVortexesParams.I2.clear();
    wakeVortexesParams.I2.resize(W.getWake().vtx.size(), { 0.0, 0.0 });
 
    wakeVortexesParams.I3.clear();
    wakeVortexesParams.I3.resize(W.getWake().vtx.size(), { 0.0, 0.0 });
 
    wakeVortexesParams.diffVelo.clear();
    wakeVortexesParams.diffVelo.resize(W.getWake().vtx.size(), { 0.0, 0.0 });
 
    wakeVortexesParams.epsastWake.clear();
    wakeVortexesParams.epsastWake.resize(W.getWake().vtx.size(), 0.0);
    W.getTimers().stop("DiffVel");
 
    //Создаем массивы под виртуальные вихри
    W.getTimers().start("ConvVel");
    virtualVortexesParams.clear();
    virtualVortexesParams.resize(W.getNumberOfBoundary());
    W.getTimers().stop("ConvVel");
 
    for (size_t bou = 0; bou < W.getNumberOfBoundary(); ++bou)
    {
        W.getTimers().start("ConvVel");
        virtualVortexesParams[bou].convVelo.clear();
        virtualVortexesParams[bou].convVelo.resize(W.getBoundary(bou).virtualWake.vtx.size(), { 0.0, 0.0 });
        W.getTimers().stop("ConvVel");
 
        W.getTimers().start("DiffVel");
        virtualVortexesParams[bou].I0.clear();
        virtualVortexesParams[bou].I0.resize(W.getBoundary(bou).virtualWake.vtx.size(), 0.0);
 
        virtualVortexesParams[bou].I1.clear();
        virtualVortexesParams[bou].I1.resize(W.getBoundary(bou).virtualWake.vtx.size(), 0.0);
 
        virtualVortexesParams[bou].I2.clear();
        virtualVortexesParams[bou].I2.resize(W.getBoundary(bou).virtualWake.vtx.size(), { 0.0, 0.0 });
 
        virtualVortexesParams[bou].I3.clear();
        virtualVortexesParams[bou].I3.resize(W.getBoundary(bou).virtualWake.vtx.size(), { 0.0, 0.0 });
 
        virtualVortexesParams[bou].diffVelo.clear();
        virtualVortexesParams[bou].diffVelo.resize(W.getBoundary(bou).virtualWake.vtx.size(), { 0.0, 0.0 });
 
        virtualVortexesParams[bou].epsastWake.clear();
        virtualVortexesParams[bou].epsastWake.resize(W.getBoundary(bou).virtualWake.vtx.size(), 0.0);
 
        W.getTimers().stop("DiffVel");
    }
}//ResizeAndZero()

Here is the call graph for this function:

SaveVisStress()

void Velocity::SaveVisStress ( )

inherited

Сохранение вязких напряжений

Вызывается в CalcDiffVelo()

Definition at line 587 of file Velocity2D.cpp.

{
    if ((W.getPassport().timeDiscretizationProperties.saveVtxStep > 0) && (W.ifDivisible(W.getPassport().timeDiscretizationProperties.saveVisStress)))
        //      if ((W.getPassport().timeDiscretizationProperties.saveVTK > 0) && (W.ifDivisible(10)) && (W.getNumberOfAirfoil() > 0))
    {
        for (size_t q = 0; q < W.getNumberOfAirfoil(); ++q)
        {
            if (q == 0)
                VMlib::CreateUserDirectory(W.getPassport().dir, "visStress");
 
            std::stringstream ss;
            ss << "VisStress_" << q << "-";
            std::string fname = VMlib::fileNameStep(ss.str(), W.getPassport().timeDiscretizationProperties.nameLength, W.getCurrentStep(), "txt");
            std::ofstream outfile;
            outfile.open(W.getPassport().dir + "visStress/" + fname);
 
            outfile << W.getAirfoil(q).viscousStress.size() << std::endl; //Сохранение числа вихрей в пелене
 
            for (size_t i = 0; i < W.getAirfoil(q).viscousStress.size(); ++i)
            {
                const Point2D& r = 0.5 * (W.getAirfoil(q).getR(i + 1) + W.getAirfoil(q).getR(i));
                double gi = W.getAirfoil(q).viscousStress[i];
                outfile << static_cast<int>(i) << " " << r[0] << " " << r[1] << " " << gi << std::endl;
            }//for i    
            outfile.close();
        }
    }
}//SaveVisStress()

Here is the call graph for this function:

Member Data Documentation

optimizedVelocity

OptimizedVelocity VM2D::Velocity::optimizedVelocity

inherited

Definition at line 117 of file Velocity2D.h.

virtualVortexesParams

std::vector<VortexesParams> VM2D::Velocity::virtualVortexesParams

inherited

Вектор струтур, определяющий параметры виртуальных вихрей для профилей

Definition at line 115 of file Velocity2D.h.

W

const World2D& VM2D::Velocity::W

protectedinherited

Константная ссылка на решаемую задачу

Definition at line 108 of file Velocity2D.h.

wakeVortexesParams

VortexesParams VM2D::Velocity::wakeVortexesParams

inherited

Струтура, определяющая параметры вихрей в следе

Definition at line 112 of file Velocity2D.h.

---

The documentation for this class was generated from the following files:

• VM2D/Velocity2D/Velocity2DBiotSavart.h
• VM2D/Velocity2D/Velocity2DBiotSavart.cpp