Commit 5f785540 authored by Tomáš Jakubec's avatar Tomáš Jakubec
Browse files

Debug.h exporting static variables from static functions

parent 64f4ce82

Unstructured_mesh/Unstructured_mesh.pro

+3 −1
Original line number Diff line number Diff line
@@ -5,6 +5,7 @@ CONFIG -= qt 


SOURCES += \

        main.cpp \

    multiphaseflow.cpp



HEADERS += \

    ../src/Macros/MacroForEach.h \
@@ -49,7 +50,8 @@ HEADERS += \ 
    ../src/UnstructuredMesh/UnstructedMeshDefine.h \

    ../src/UnstructuredMesh/UnstructuredMesh.h \

    ../src/NumericStaticArray/Vector.h \

    ../src/NumericStaticArray/Vertex.h

    ../src/NumericStaticArray/Vertex.h \

    multiphaseflow.h



DISTFILES += \

    ../src/Debug/README.md \

Unstructured_mesh/main.cpp

+130 −213
Original line number Diff line number Diff line
@@ -2,144 +2,138 @@ 
//#define UNDEBUG

#define CONSOLE_COLORED_OUTPUT

#include "../src/Debug/Debug.h"

#include "../src/UnstructuredMesh/UnstructuredMesh.h"

#include "../src/UnstructuredMesh/MeshFunctions/MeshFunctions.h"

#include "../src/UnstructuredMesh/MeshIO/MeshReader/VTKMeshReader.h"

#include "../src/UnstructuredMesh/MeshIO/MeshWriter/VTKMeshWriter.h"

#include "../src/UnstructuredMesh/MeshDataContainer/MeshDataIO/VTKMeshDataWriter.h"

#include "../src/UnstructuredMesh/MeshDataContainer/MeshDataIO/VTKMeshDataReader.h"

#include "../src/UnstructuredMesh/MeshIO/MeshReader/FPMAMeshReader.h"

#include "../src/UnstructuredMesh/MeshIO/MeshWriter/FPMAMeshWriter.h"



#include "../src/Traits/Traits.h"

#include "../src/Traits/TraitsAlgorithm/TraitsAlgorithm.h"

#include "multiphaseflow.h"



#include <fstream>

#include <list>

#include <chrono>

using namespace std;







struct CellData {

    double invVol;



};





struct CompData {

    double T;

    //CompData(const CompData&) = default;

    //CompData(CompData&&) = default;

    CompData(double _T = 0){T = _T;}



};

MAKE_NAMED_ATTRIBUTE_TRAIT(CompData, "Temperature", T);

MAKE_ATTRIBUTE_TRAIT_ARITHMETIC(CompData, T);



struct FaceData {

    double volOverDist;

};



template <unsigned int ProblemDimension, typename Real>

class HeatCunduction {

public:

    using MeshType = UnstructuredMesh<ProblemDimension,size_t, double, 12>;

    using ResultType = CompData;



public:

    MeshType mesh;

    MeshDataContainer<std::tuple<CellData, FaceData>, 3,2> meshData;

    VTKMeshWriter<ProblemDimension, size_t, Real> writer;



    HeatCunduction(){

    }



    void loadMesh(const std::string& fileName) {

        FPMAMeshReader<MeshType::meshDimension()> reader;

        ifstream file(fileName);

        reader.loadFromStream(file, mesh);



        mesh.template initializeCenters<TESSELLATED>();



        mesh.setupBoundaryCells();







        mesh.setupBoundaryCellsCenters();





        auto measures = mesh.template computeElementMeasures<TESSELLATED>();



        auto dists = ComputeCellsDistance(mesh);





        meshData.allocateData(mesh);



        MeshDataContainer<CompData, 3> compData(mesh);



        for (auto& face : mesh.getFaces()) {

            meshData.at(face).volOverDist = measures.at(face) / dists.at(face);

        }



        for (auto& cell : mesh.getCells()) {

            meshData.at(cell).invVol = 1.0 / measures.at(cell);

        }

    }



    void calculateRHS(

            Real,//time is unused in this problem

            MeshDataContainer<CompData, ProblemDimension>& compData,

            MeshDataContainer<CompData, ProblemDimension>& outDeltas) {



        for (CompData& dT : outDeltas.template getDataByPos<0>()){

            dT.T = 0;

        }



        for (auto& face : mesh.getFaces()) {

            CompData lT(1000);

            bool lIn = false;

            if(!((face.getCellLeftIndex() & BOUNDARY_INDEX(size_t)) == BOUNDARY_INDEX(size_t))){

                lT = compData.template getDataByDim<3>().at(face.getCellLeftIndex());

                lIn = true;

            }



            CompData rT(1000);

            bool rIn = false;

            if(!((face.getCellRightIndex() & BOUNDARY_INDEX(size_t)) == BOUNDARY_INDEX(size_t))){

                rT = compData.template getDataByDim<3>().at(face.getCellRightIndex());

                rIn = true;

            }





            CompData dT = meshData.at(face).volOverDist * (lT - rT);

            if (rIn)

                outDeltas.template getDataByDim<3>().at(face.getCellRightIndex()) += dT;

            if (lIn)

                outDeltas.template getDataByDim<3>().at(face.getCellLeftIndex()) -= dT;



        }



        for (auto& cell : mesh.getCells()) {

            outDeltas.at(cell) *= meshData.at(cell).invVol;

        }

    }



    template<typename dataType>

    void exportData(double time,

                    MeshDataContainer<dataType, 3>& compData) {

using namespace std;



        ofstream ofile("HeatCond"s + "_" + to_string(time) + ".vtk");

        writer.writeHeader(ofile, "HC_test"s + to_string(time));

        writer.writeToStream(ofile, mesh, MeshDataContainer<MeshNativeType<3>::ElementType, 3>(mesh, MeshNativeType<3>::POLYHEDRON));



        VTKMeshDataWriter<3> dataWriter;

        dataWriter.writeToStream(ofile, compData, writer);



        ofile.close();

        DBGMSG("Data eported");



    }

};

// struct CellData {

//     double invVol;

//

// };

//

//

// struct CompData {

//     double T;

//     //CompData(const CompData&) = default;

//     //CompData(CompData&&) = default;

//     CompData(double _T = 0){T = _T;}

//

// };

// MAKE_NAMED_ATTRIBUTE_TRAIT(CompData, "Temperature", T);

// MAKE_ATTRIBUTE_TRAIT_ARITHMETIC(CompData, T);

//

// struct FaceData {

//     double volOverDist;

// };

//

// template <unsigned int ProblemDimension, typename Real>

// class HeatCunduction {

// public:

//     using MeshType = UnstructuredMesh<ProblemDimension,size_t, double, 12>;

//     using ResultType = CompData;

//

// public:

//     MeshType mesh;

//     MeshDataContainer<std::tuple<CellData, FaceData>, 3,2> meshData;

//     VTKMeshWriter<ProblemDimension, size_t, Real> writer;

//

//     HeatCunduction(){

//     }

//

//     void loadMesh(const std::string& fileName) {

//         FPMAMeshReader<MeshType::meshDimension()> reader;

//         ifstream file(fileName);

//         reader.loadFromStream(file, mesh);

//

//         mesh.template initializeCenters<TESSELLATED>();

//

//         mesh.setupBoundaryCells();

//

//

//

//         mesh.setupBoundaryCellsCenters();

//

//

//         auto measures = mesh.template computeElementMeasures<TESSELLATED>();

//

//         auto dists = ComputeCellsDistance(mesh);

//

//

//         meshData.allocateData(mesh);

//

//         MeshDataContainer<CompData, 3> compData(mesh);

//

//         for (auto& face : mesh.getFaces()) {

//             meshData.at(face).volOverDist = measures.at(face) / dists.at(face);

//         }

//

//         for (auto& cell : mesh.getCells()) {

//             meshData.at(cell).invVol = 1.0 / measures.at(cell);

//         }

//     }

//

//     void calculateRHS(

//             Real,//time is unused in this problem

//             MeshDataContainer<CompData, ProblemDimension>& compData,

//             MeshDataContainer<CompData, ProblemDimension>& outDeltas) {

//

//         for (CompData& dT : outDeltas.template getDataByPos<0>()){

//             dT.T = 0;

//         }

//

//         for (auto& face : mesh.getFaces()) {

//             CompData lT(1000);

//             bool lIn = false;

//             if(!((face.getCellLeftIndex() & BOUNDARY_INDEX(size_t)) == BOUNDARY_INDEX(size_t))){

//                 lT = compData.template getDataByDim<3>().at(face.getCellLeftIndex());

//                 lIn = true;

//             }

//

//             CompData rT(1000);

//             bool rIn = false;

//             if(!((face.getCellRightIndex() & BOUNDARY_INDEX(size_t)) == BOUNDARY_INDEX(size_t))){

//                 rT = compData.template getDataByDim<3>().at(face.getCellRightIndex());

//                 rIn = true;

//             }

//

//

//             CompData dT = meshData.at(face).volOverDist * (lT - rT);

//             if (rIn)

//                 outDeltas.template getDataByDim<3>().at(face.getCellRightIndex()) += dT;

//             if (lIn)

//                 outDeltas.template getDataByDim<3>().at(face.getCellLeftIndex()) -= dT;

//

//         }

//

//         for (auto& cell : mesh.getCells()) {

//             outDeltas.at(cell) *= meshData.at(cell).invVol;

//         }

//     }

//

//     template<typename dataType>

//     void exportData(double time,

//                     MeshDataContainer<dataType, 3>& compData) {

//

//         using std::string_literals;

//

//         ofstream ofile("HeatCond"s + "_" + to_string(time) + ".vtk");

//         writer.writeHeader(ofile, "HC_test"s + to_string(time));

//         writer.writeToStream(ofile, mesh, MeshDataContainer<MeshNativeType<3>::ElementType, 3>(mesh, MeshNativeType<3>::POLYHEDRON));

//

//         VTKMeshDataWriter<3> dataWriter;

//         dataWriter.writeToStream(ofile, compData, writer);

//

//         ofile.close();

//         DBGMSG("Data eported");

//

//     }

// };



struct watch{
@@ -284,96 +278,19 @@ wError.lap(); 






struct Index {

    size_t index;

};



MAKE_ATTRIBUTE_TRAIT(Index, index);



void meshAdmisibility(const std::string& meshFileName) {

    UnstructuredMesh<3, size_t, double, 12> mesh;

    FPMAMeshReader<3> reader;

    ifstream file(meshFileName);

    reader.loadFromStream(file, mesh);



    mesh.initializeCenters<TESSELLATED>();



    mesh.setupBoundaryCells();







    mesh.setupBoundaryCellsCenters();





    auto faceNormals = mesh.computeFaceNormals<TESSELLATED>();



    MeshDataContainer<Vector<3, double>, 2> centerLines(mesh);



    double max = -1, min = 2, avg = 0;

    MeshDataContainer<double, 2> faceAdmisibility(mesh);



    auto faceVert = MeshConnections<2,0, ORDER_ORIGINAL>::connections(mesh);



    int cnt = 0;

    for( auto& face : mesh.getFaces()) {

        //DBGVAR(faceIndex);

        bool badFace = false;

        for (size_t i = 1; i < faceVert.at(face).size(); i++) {

            size_t index = faceVert.at(face).at(i);

            Vertex<3,double> centerToVert = mesh.getVertices().at(faceVert.at(face).at(0)) - mesh.getVertices().at(index);

            double scalarProduct = faceNormals.at(face) * (Vector<3,double>{centerToVert[0],centerToVert[1],centerToVert[2]});

            //HTMLDBGCOND(scalarProduct > 0.01, faceNormals.at(face), centerToVert, scalarProduct);

            if (abs(scalarProduct) > 1e-4) badFace = true; // the face is not planar

        }

        if (badFace) cnt++;

    }

    DBGVAR(cnt, mesh.getFaces().size());



    for (auto& face : mesh.getFaces()) {



        auto& cellLeft = (face.getCellLeftIndex() & BOUNDARY_INDEX(size_t)) == BOUNDARY_INDEX(size_t)?

                  mesh.getBoundaryCells().at(face.getCellLeftIndex()&EXTRACTING_INDEX(size_t)):

                  mesh.getCells().at(face.getCellLeftIndex());

        auto& cellRight = (face.getCellRightIndex() & BOUNDARY_INDEX(size_t)) == BOUNDARY_INDEX(size_t)?

                  mesh.getBoundaryCells().at(face.getCellRightIndex()&EXTRACTING_INDEX(size_t)):

                  mesh.getCells().at(face.getCellRightIndex());





        Vertex<3,double> line = cellRight.getCenter() - cellLeft.getCenter();

        centerLines.at(face) = {line[0],line[1],line[2]};

        centerLines.at(face) /= centerLines.at(face).normEukleid();



        faceAdmisibility.at(face) = faceNormals.at(face) * centerLines.at(face);

        if (faceAdmisibility.at(face) > max) {

            max = faceAdmisibility.at(face);

        }

        if (faceAdmisibility.at(face) < min) {

            min = faceAdmisibility.at(face);

        }

        avg += faceAdmisibility.at(face);

        //DBGVARCOND((face.getCellLeftIndex() == 691 || face.getCellLeftIndex() == 699) && (face.getCellRightIndex() == 691 || face.getCellRightIndex() == 699), face.getIndex());

        //DBGVARCOND_HTML(faceAdmisibility.at(face) < 0.0, faceNormals.at(face), centerLines.at(face), face.getIndex(), face.getCenter(), face.getCellLeftIndex(), cellLeft.getCenter(), face.getCellRightIndex(), cellRight.getCenter(), faceAdmisibility.at(face));

    }

    avg /= mesh.getFaces().size();

    DBGVAR(max, min, avg);



    VTKMeshWriter<3, size_t, double> writer;

    MeshDataContainer<Index, 3> indexes(mesh);

    for (size_t i = 0; i < mesh.getCells().size(); i++) {

        indexes.getDataByDim<3>().at(i).index = i;

    }

    //exportData("Admisibility", writer, 0, mesh, indexes);

}









void testHeatConduction1() {



    MultiphaseFlow mpf;



void testHeatConduction1() {

    mpf.setupMeshData("Boiler.vtk");



/*

    HeatCunduction<3,double> hcProblem;



    hcProblem.loadMesh("Poly_2_5_level2.fpma");



    double startTime = 0.0, finalTime = 0.1 , tau = 1e-4;
@@ -408,7 +325,7 @@ wK1.reset();wK2.reset();wK3.reset();wK4.reset();wError.reset(); 
//    hcProblem.exportData(finalTime, compData);



DBGVAR(wK1.getResult(),wK2.getResult(),wK3.getResult(),wK4.getResult(),wError.getResult());



*/





}

Unstructured_mesh/multiphaseflow.cpp

0 → 100644
+503 −0

File added.

Preview size limit exceeded, changes collapsed.

Unstructured_mesh/multiphaseflow.h

0 → 100644
+913 −0

File added.

Preview size limit exceeded, changes collapsed.

src/Debug/Debug.h

+12 −7
Original line number Diff line number Diff line
@@ -14,12 +14,17 @@ 
*/

namespace dbg {

    struct DBGStatics {

        static HtmlLogger HDBGLog;

        static CSVLogger CSVDBGLog;

    };

        static HtmlLogger getHTMLLogger(){

            static HtmlLogger HDBGLog("DBG.html");

            return HDBGLog;

        }





    HtmlLogger DBGStatics::HDBGLog("DBG.html");

    CSVLogger DBGStatics::CSVDBGLog("DBG.csv");

        static CSVLogger getCSVLogger(){

            static CSVLogger CSVDBGLog("DBG.csv");

            return CSVDBGLog;

        }

    };

}



#define STRVAR(var) #var, var
@@ -37,11 +42,11 @@ ConsoleLogger<>::writeMessage("!!", __LINE__, __FILE__, std::string("something w 
abort();}



// Macros using html debug output

#define DBGVAR_HTML(...) dbg::DBGStatics::HDBGLog.writeVar(__LINE__, __FILE__, FOR_EACH(STRVAR, __VA_ARGS__))

#define DBGVAR_HTML(...) dbg::DBGStatics::getHTMLLogger().writeVar(__LINE__, __FILE__, FOR_EACH(STRVAR, __VA_ARGS__))

#define DBGVARCOND_HTML(condition, ...) if(condition) DBGVAR_HTML(__VA_ARGS__)



// Macros using csv debug output

#define DBGVAR_CSV(...) dbg::DBGStatics::CSVDBGLog.writeVar(__LINE__, __FILE__, FOR_EACH(STRVAR, __VA_ARGS__))

#define DBGVAR_CSV(...) dbg::DBGStatics::getCSVLogger().writeVar(__LINE__, __FILE__, FOR_EACH(STRVAR, __VA_ARGS__))

#define DBGVARCOND_CSV(condition, ...) if(condition) DBGVAR_HTML(__VA_ARGS__)