PCM.h

/*-------------------------------------------------------------------
Copyright 2013 Ravishankar Sundararaman, Deniz Gunceler

This file is part of JDFTx.

JDFTx 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, either version 3 of the License, or
(at your option) any later version.

JDFTx 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.

You should have received a copy of the GNU General Public License
along with JDFTx.  If not, see <http://www.gnu.org/licenses/>.
-------------------------------------------------------------------*/

#ifndef JDFTX_ELECTRONIC_PCM_H
#define JDFTX_ELECTRONIC_PCM_H

#include <fluid/FluidSolver.h>
#include <electronic/RadialFunction.h>
#include <core/EnergyComponents.h>
#include <core/Coulomb.h>

class PCM : public FluidSolver
{
public:
        PCM(const Everything& e, const FluidSolverParams& fsp);
        virtual ~PCM();
        
        void dumpDensities(const char* filenamePattern) const; 
        void dumpDebug(const char* filenamePattern) const; 

protected:
        EnergyComponents Adiel; 
        ScalarFieldTilde rhoExplicitTilde; 
        ScalarField nCavity, tauCavity, nCavityEx[2]; 
        ScalarField shape, shapeVdw; 
        
        virtual void printDebug(FILE* fp) const {} 
        
        void updateCavity(); 
        void propagateCavityGradients(const ScalarField& A_shape, ScalarField& A_nCavity, ScalarFieldTilde& A_rhoExplicitTilde, bool electricOnly) const; 
        void setExtraForces(IonicGradient* forces, const ScalarFieldTilde& A_nCavityTilde) const; 
        ScalarFieldTilde getFullCore() const; 
private:
        ScalarField Acavity_shape, Acavity_shapeVdw; 
        double A_nc, A_tension, A_vdwScale, A_eta_wDiel, A_pCavity; 
        double Rex[2]; 
        RadialFunctionG wExpand[2]; 
        RadialFunctionG wCavity; 
protected:
        std::vector<RadialFunctionG> Sf; 
        std::vector<int> atomicNumbers; 
        static ScalarFieldTilde coulomb(const ScalarFieldTilde& rho) { return (-4*M_PI) * Linv(O(rho)); }
        friend struct ChargedDefect;
};

#endif // JDFTX_ELECTRONIC_PCM_H