PFSolver.h

/* Copyright 2017-2021 Institute for Automation of Complex Power Systems,
 *                     EONERC, RWTH Aachen University
 *
 * This Source Code Form is subject to the terms of the Mozilla Public
 * License, v. 2.0. If a copy of the MPL was not distributed with this
 * file, You can obtain one at https://mozilla.org/MPL/2.0/.
 *********************************************************************************/
 
#pragma once
 
#include <cmath>
#include <iterator>
 
#include "dpsim-models/Components.h"
#include "dpsim-models/SystemTopology.h"
#include <dpsim/Scheduler.h>
#include <dpsim/Solver.h>
 
namespace DPsim {
class PFSolver : public Solver {
protected:
  UInt mNumPQBuses = 0;
  UInt mNumPVBuses = 0;
  UInt mNumVDBuses = 0;
  UInt mNumUnknowns = 0;
  CPS::TopologicalNode::List mPQBuses;
  CPS::TopologicalNode::List mPVBuses;
  CPS::TopologicalNode::List mVDBuses;
  CPS::TopologicalNode::List mPQBusesOrig;
  CPS::TopologicalNode::List mPVBusesOrig;
  std::vector<CPS::UInt> mPQBusIndices;
  std::vector<CPS::UInt> mPVBusIndices;
  std::vector<CPS::UInt> mVDBusIndices;
  std::vector<CPS::UInt> mPQPVBusIndices;

  CPS::SparseMatrixCompRow mY;

  CPS::Matrix mJ;
  CPS::Vector mX;
  CPS::Vector mF;

  CPS::SystemTopology mSystem;
  std::vector<std::shared_ptr<CPS::SP::Ph1::Transformer>> mTransformers;
  std::vector<std::shared_ptr<CPS::SP::Ph1::SolidStateTransformer>>
      mSolidStateTransformers;
  std::vector<std::shared_ptr<CPS::SP::Ph1::SynchronGenerator>>
      mSynchronGenerators;
  std::vector<std::shared_ptr<CPS::SP::Ph1::Load>> mLoads;
  std::vector<std::shared_ptr<CPS::SP::Ph1::PiLine>> mLines;
  std::vector<std::shared_ptr<CPS::SP::Ph1::Shunt>> mShunts;
  std::vector<std::shared_ptr<CPS::SP::Ph1::NetworkInjection>> mExternalGrids;
  std::vector<std::shared_ptr<CPS::SP::Ph1::AvVoltageSourceInverterDQ>>
      mAverageVoltageSourceInverters;
  std::map<CPS::TopologicalNode::Ptr, CPS::Real> mBaseVoltageAtNode;

  Real mTolerance = 1e-8;
  CPS::UInt mMaxIterations = 20;
  CPS::UInt mIterations;
  CPS::Bool mEnforceReactiveLimits = false;
  CPS::UInt mMaxOuterIterations = 10;
  CPS::UInt mMaxQLimitSwitchesPerBus = 2;
  CPS::Real mBaseVoltageLooseTolerance = 0.1;
  CPS::Real mBaseVoltageStrictTolerance = 0.01;
  CPS::Real mBaseApparentPower;
  CPS::Real mBaseApparentPowerFallback = 100e6;
  CPS::Bool isConverged = false;
  CPS::Bool solutionInitialized = false;
  CPS::Bool solutionComplexInitialized = false;

  CPS::Bool mKeepLastSolution = false;

  virtual void generateInitialSolution(Real time,
                                       bool keep_last_solution = false) = 0;
  virtual void calculateMismatch() = 0;
  virtual void calculateJacobian() = 0;
  virtual void updateSolution() = 0;
  virtual void setSolution() = 0;

  void initialize() override;
  void initializeComponents();
  void assignMatrixNodeIndices();
  void setBaseApparentPower();
  void determinePFBusType();
  void rebuildBusIndexAggregates();
  void reclassifyBuses();
  void resetToOriginalClassification();
  virtual void clearReactiveLimitState() {}
  CPS::Real componentBaseVoltage(CPS::TopologicalPowerComp::Ptr comp,
                                 CPS::TopologicalNode::Ptr node);
  void propagateAndVerifyBaseVoltage();

  void composeAdmittanceMatrix();
  CPS::Real G(int i, int j);
  CPS::Real B(int i, int j);
  Bool solvePowerflow();
  Bool runNewtonRaphson();
  virtual CPS::Bool enforceReactiveLimits() { return false; }
  virtual void setUpJacobianStorage();
  virtual void solveJacobianSystem();
  CPS::Bool checkConvergence();
  CPS::String logVector(std::vector<CPS::UInt> indexVector) {
    std::stringstream result;
    std::copy(indexVector.begin(), indexVector.end(),
              std::ostream_iterator<CPS::UInt>(result, " "));
    return result.str();
  };

  CPS::Task::List getTasks() override;
  // determines power flow bus type for each node according to the components attached to it.
public:
  PFSolver(CPS::String name, CPS::SystemTopology system, Real timeStep,
           CPS::Logger::Level logLevel);
  virtual ~PFSolver(){};

  void setVDNode(CPS::String name);
  void modifyPowerFlowBusComponent(CPS::String name,
                                   CPS::PowerflowBusType powerFlowBusType);
  void setSolverAndComponentBehaviour(Solver::Behaviour behaviour) override;
 
  void setKeepLastSolution(CPS::Bool keepLastSolution) {
    mKeepLastSolution = keepLastSolution;
  }

  void setEnforceReactiveLimits(CPS::Bool value) {
    mEnforceReactiveLimits = value;
  }

  void setBaseVoltageLooseTolerance(CPS::Real tolerance) {
    mBaseVoltageLooseTolerance = tolerance;
  }

  void setBaseVoltageStrictTolerance(CPS::Real tolerance) {
    mBaseVoltageStrictTolerance = tolerance;
  }
 
  CPS::Bool getKeepLastSolution() const { return mKeepLastSolution; }
 
  void setBaseApparentPowerFallback(CPS::Real baseApparentPowerFallback) {
    mBaseApparentPowerFallback = baseApparentPowerFallback;
  }
 
  CPS::Real getBaseApparentPowerFallback() const {
    return mBaseApparentPowerFallback;
  }
 
  void setMaxIterations(CPS::UInt maxIterations) {
    mMaxIterations = maxIterations;
  }
 
  CPS::UInt getMaxIterations() const { return mMaxIterations; }
 
  class SolveTask : public CPS::Task {
  public:
    SolveTask(PFSolver &solver)
        : Task(solver.mName + ".Solve"), mSolver(solver) {
      mModifiedAttributes.push_back(Scheduler::external);
    }
 
    void execute(Real time, Int timeStepCount);
 
  private:
    PFSolver &mSolver;
  };
};
} // namespace DPsim