electr2d.hpp

Go to the documentation of this file.

/*
 * This file is part of PLaSK (https://plask.app) by Photonics Group at TUL
 * Copyright (c) 2022 Lodz University of Technology
 *
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, version 3.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
 * GNU General Public License for more details.
 */
#ifndef PLASK__MODULE_ELECTRICAL_ELECTR2D_H
#define PLASK__MODULE_ELECTRICAL_ELECTR2D_H
 
#include "common.hpp"
 
namespace plask { namespace electrical { namespace shockley {
 
template <typename Geometry2DType>
struct PLASK_SOLVER_API ElectricalFem2DSolver : public FemSolverWithMaskedMesh<Geometry2DType, RectangularMesh<2>> {
  protected:
 
    struct Active {
        struct Region {
            size_t left, right, bottom, top;
            size_t rowl, rowr;
            bool warn;
            Region()
                : left(0),
                  right(0),
                  bottom(std::numeric_limits<size_t>::max()),
                  top(std::numeric_limits<size_t>::max()),
                  rowl(std::numeric_limits<size_t>::max()),
                  rowr(0),
                  warn(true) {}
        };
        size_t left, right, bottom, top;
        ptrdiff_t offset;
        double height;
        Active() : left(0), right(0), bottom(0), top(0), offset(0), height(0.) {}
        Active(size_t tot, size_t l, size_t r, size_t b, size_t t, double h)
            : left(l), right(r), bottom(b), top(t), offset(tot - l), height(h) {}
    };
 
    double pcond;  
    double ncond;  
 
    int loopno;     
    double toterr;  
    Vec<2, double> maxcur;  
 
    DataVector<Tensor2<double>> junction_conductivity;  
    Tensor2<double> default_junction_conductivity;      
 
    DataVector<Tensor2<double>> conds;    
    DataVector<double> potentials;        
    DataVector<Vec<2, double>> currents;  
    DataVector<double> heats;             
 
    std::vector<Active> active;  
 
    inline void setLocalMatrix(double& k44,
                               double& k33,
                               double& k22,
                               double& k11,
                               double& k43,
                               double& k21,
                               double& k42,
                               double& k31,
                               double& k32,
                               double& k41,
                               double ky,
                               double width,
                               const Vec<2, double>& midpoint);
 
    virtual Tensor2<double> activeCond(size_t n, double U, double jy, double T) = 0;
 
    LazyData<double> loadConductivities();
 
    void saveConductivities();
 
    void saveHeatDensities();
 
    void onInitialize() override;
 
    void onInvalidate() override;
 
    void setupActiveRegions();
 
    // void onMeshChange(const typename RectangularMesh<2>::Event& evt) override {
    //     SolverWithMesh<Geometry2DType, RectangularMesh<2>>::onMeshChange(evt);
    //     setupActiveRegions();
    // }
 
    // void onGeometryChange(const Geometry::Event& evt) override {
    //     SolverWithMesh<Geometry2DType, RectangularMesh<2>>::onGeometryChange(evt);
    //     setupActiveRegions();
    // }
 
    void setMatrix(FemMatrix& A,
                   DataVector<double>& B,
                   const BoundaryConditionsWithMesh<RectangularMesh<2>::Boundary, double>& bvoltage,
                   const LazyData<double>& temperature);
 
    size_t isActive(const Vec<2>& point) const {
        size_t no(0);
        auto roles = this->geometry->getRolesAt(point);
        for (auto role : roles) {
            size_t l = 0;
            if (role.substr(0, 6) == "active")
                l = 6;
            else if (role.substr(0, 8) == "junction")
                l = 8;
            else
                continue;
            if (no != 0) throw BadInput(this->getId(), "multiple 'active'/'junction' roles specified");
            if (role.size() == l)
                no = 1;
            else {
                try {
                    no = boost::lexical_cast<size_t>(role.substr(l)) + 1;
                } catch (boost::bad_lexical_cast&) {
                    throw BadInput(this->getId(), "bad junction number in role '{0}'", role);
                }
            }
        }
        return no;
    }
 
    size_t isActive(const RectangularMaskedMesh2D::Element& element) const { return isActive(element.getMidpoint()); }
 
  public:
    double maxerr;  
 
    BoundaryConditions<RectangularMesh<2>::Boundary, double> voltage_boundary;
 
    typename ProviderFor<Voltage, Geometry2DType>::Delegate outVoltage;
 
    typename ProviderFor<CurrentDensity, Geometry2DType>::Delegate outCurrentDensity;
 
    typename ProviderFor<Heat, Geometry2DType>::Delegate outHeat;
 
    typename ProviderFor<Conductivity, Geometry2DType>::Delegate outConductivity;
 
    ReceiverFor<Temperature, Geometry2DType> inTemperature;
 
    Convergence convergence;    
 
    double compute(unsigned loops = 1);
 
    double integrateCurrent(size_t vindex, bool onlyactive = false);
 
    double getTotalCurrent(size_t nact = 0);
 
    double getTotalEnergy();
 
    double getCapacitance();
 
    double getTotalHeat();
 
    double getErr() const { return toterr; }
 
    double getCondPcontact() const { return pcond; }
    void setCondPcontact(double cond) {
        pcond = cond;
        this->invalidate();
    }
 
    double getCondNcontact() const { return ncond; }
    void setCondNcontact(double cond) {
        ncond = cond;
        this->invalidate();
    }
 
    DataVector<const Tensor2<double>> getCondJunc() const { return junction_conductivity; }
    void setCondJunc(double cond) {
        junction_conductivity.reset(max(junction_conductivity.size(), size_t(1)), cond);
        default_junction_conductivity = Tensor2<double>(0., cond);
    }
    void setCondJunc(Tensor2<double> cond) {
        junction_conductivity.reset(max(junction_conductivity.size(), size_t(1)), cond);
        default_junction_conductivity = cond;
    }
    void setCondJunc(const DataVector<Tensor2<double>>& cond) {
        size_t condsize = 0;
        condsize = max(condsize, size_t(1));
        if (!this->mesh || cond.size() != condsize)
            throw BadInput(this->getId(), "provided junction conductivity vector has wrong size");
        junction_conductivity = cond.claim();
    }
 
    void loadConfiguration(XMLReader& source, Manager& manager) override;
 
    void parseConfiguration(XMLReader& source, Manager& manager);
 
    ElectricalFem2DSolver(const std::string& name = "");
 
    ~ElectricalFem2DSolver();
 
  protected:
    const LazyData<double> getVoltage(shared_ptr<const MeshD<2>> dest_mesh, InterpolationMethod method) const;
 
    const LazyData<double> getHeatDensities(shared_ptr<const MeshD<2>> dest_mesh, InterpolationMethod method);
 
    const LazyData<Vec<2>> getCurrentDensities(shared_ptr<const MeshD<2>> dest_mesh, InterpolationMethod method);
 
    const LazyData<Tensor2<double>> getConductivity(shared_ptr<const MeshD<2>> dest_mesh, InterpolationMethod method);
};
 
}}}  // namespace plask::electrical::shockley
 
#endif