diff --git a/docs/developer-guide/andes-compatibility.md b/docs/developer-guide/andes-compatibility.md new file mode 100644 index 000000000..72d5d6eb2 --- /dev/null +++ b/docs/developer-guide/andes-compatibility.md @@ -0,0 +1,98 @@ +# ANDES compatibility roadmap + +This document records the work required for GridDyn to load and solve ANDES +cases. It is intentionally a living checklist: a model is not complete until +its importer mapping and a numerical regression test are both present. + +## Status terms + +| Status | Meaning | +| ------------------ | -------------------------------------------------------------------------------------------- | +| Implemented | GridDyn model, ANDES importer mapping, and a numerical test exist. | +| Partial | A GridDyn analogue exists, but the mapping, controls, or numerical validation is incomplete. | +| Planned | No implementation work has started. | +| Untriaged | Included in the ANDES inventory, but its GridDyn-equivalence review has not started. | +| No direct analogue | Requires a new GridDyn model rather than reader-only work. | + +## Power-flow model mapping + +| ANDES model(s) | GridDyn mapping | Status | Notes | +| -------------------------------------------------- | -------------------------- | ------------------ | ------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------------- | +| `Bus` | `AcBus` | Implemented | Base voltage, initial voltage, and angle are imported. | +| `PQ` | `ZipLoad` | Implemented | Constant-power portion is imported. | +| `PV`, `Slack` | `AcBus` plus `Generator` | Implemented | Active-power and voltage targets are imported. | +| `Line` | `AcLine` | Partial | `r`, `x`, `b`, tap, phase shift, and status are imported. ANDES line-rated-voltage (`Vn1`/`Vn2`) base conversion, including its implicit 110 kV default, still needs an equivalent mixed-base implementation. | +| `Area` | `GridArea` | Partial | A GridDyn analogue exists; ANDES `Area` is not yet imported. | +| `Shunt` | fixed-admittance `ZipLoad` | Partial | The equivalent exists but needs a direct mapping and test. | +| `ShuntSw` | `loads::Svd` | Partial | Both provide switched reactive support; block/control semantics need mapping and numerical tests. | +| `ShuntTD` | fixed shunt for power flow | Partial | Its steady-state behavior follows `Shunt`; its time-domain phase-voltage outputs are not represented. | +| `Jumper` | `links::ZBreaker` | Partial | Both merge buses through a zero-impedance connection; no importer mapping/test yet. | +| `Motor3`, `Motor5` | `MotorLoad3`, `MotorLoad5` | Partial | Candidate GridDyn models exist; parameter mapping and initialization comparisons remain. | +| `Fortescue` | none | No direct analogue | Requires a multi-terminal positive-/negative-/zero-sequence interface model. | +| `Node`, `Ground` | `DcBus` | Implemented | Ground is imported as a DC swing reference. | +| `R`, `L`, `C`, `RLs`, `RCp`, `RLCp`, `RCs`, `RLCs` | `links::DcLink` | Implemented | ANDES current-balance convention is retained. Dynamic-form validation remains open. | +| `VSCShunt` | `links::VSCShunt` | Implemented | Three-terminal algebraic converter; impedance-base conversion is applied during import. | + +## Dynamic-model roadmap + +The present numerical suite covers power flow. Dynamic compatibility requires +matching model states, initialization, events, limits, and trajectories; not +only the final power-flow point. + +| Work item | Status | Completion evidence | +| ---------------------------------------------- | ------- | ------------------------------------------------------------------------------------------------------------- | +| DC `L`/`C`/combined branch dynamics | Planned | ANDES and GridDyn trajectories from the same initialized DC case. | +| `VSCShunt` controls, losses, limits, and droop | Partial | Core algebraic controls are present; loss coefficients and limit enforcement need model and trajectory tests. | +| Motors (`Motor3`/`Motor5`) | Planned | State-name/initialization mapping plus disturbance trajectory comparisons. | +| Switched shunts | Planned | Block selection and switching-event comparisons. | +| Fortescue interface | Planned | New GridDyn component and unbalanced/interface regression cases. | +| Remaining ANDES dynamic families | Planned | Maintain a model-by-model inventory before adding importer mappings. | + +## Dynamic-model inventory + +This inventory is based on the 100 registered models in the local ANDES +installation. Rows intentionally retain the ANDES class names so that a row +can be split and checked off as soon as models in the same family diverge. A +model is not compatible merely because a similarly named GridDyn model exists: +it also needs native-input mapping, initialization, and a trajectory test. + +| ANDES model(s) | GridDyn mapping / next action | Status | +| ------------------------------------------------------------------------------------------------------------------------------------------------------ | ----------------------------------------------------------------------------------------------------------------- | --------- | +| `GENCLS` | `GenModelClassical`; parameter and trajectory comparison required. | Partial | +| `GENROU` | `GenModelGENROU`; PSS/e DYR reader has a mapping, but native ANDES mapping and numerical validation are absent. | Partial | +| `TGOV1` | GridDyn governor; PSS/e DYR reader has a mapping, but native ANDES mapping and numerical validation are absent. | Partial | +| `EXDC2` | GridDyn DC exciter; PSS/e DYR reader has a mapping, but native ANDES mapping and numerical validation are absent. | Partial | +| `ZIP`, `FLoad` | GridDyn static/dynamic load models; identify exact parameter and frequency-response equivalence. | Partial | +| `Motor3`, `Motor5` | `MotorLoad3`, `MotorLoad5`; parameter mapping and trajectory comparisons required. | Partial | +| `ACE`, `ACEc`, `COI` | Map area-control and center-of-inertia services. | Untriaged | +| `PLBVFU1`, `IEEEVC` | Map playback and voltage-compensation controls. | Untriaged | +| `TG2`, `TGOV1DB`, `TGOV1N`, `TGOV1NDB`, `IEEEG1`, `IEESGO`, `GAST`, `HYGOV`, `HYGOVDB`, `HYGOV4` | Governor-family inventory and model-by-model equivalence review. | Untriaged | +| `SHAFT5` | Multi-mass shaft model; map states and mechanical interfaces. | Untriaged | +| `IEEEX1`, `ESDC1A`, `ESDC2A`, `EXST1`, `ESST3A`, `SEXS`, `IEEET1`, `EXAC1`, `EXAC2`, `EXAC4`, `ESST4B`, `AC8B`, `IEEET3`, `ESAC1A`, `ESST1A`, `ESAC5A` | Exciter-family inventory; review the remaining classes individually. | Untriaged | +| `IEEEST`, `ST2CUT` | Power-system stabilizer models; identify GridDyn analogue or implement. | Untriaged | +| `BusFreq`, `BusROCOF`, `PMU`, `PLL1`, `PLL2`, `FreqDiv` | Measurement and frequency-estimation models. | Untriaged | +| `REGCA1`, `REGCP1`, `REECA1`, `REECA1E`, `REECA1G`, `REECB1`, `REPCA1`, `REGCV1`, `REGCV2`, `REGF1`, `REGF2`, `REGF3` | Renewable/inverter-controller family; map plant, electrical, and frequency-control blocks. | Untriaged | +| `WTDTA1`, `WTDS`, `WTARA1`, `WTPTA1`, `WTTQA1`, `WTARV1` | Wind-turbine drive-train and control models. | Untriaged | +| `PVD1`, `ESD1`, `EV1`, `EV2`, `DGPRCT1`, `DGPRCTExt` | Distributed energy-resource and protection models. | Untriaged | +| `Fault`, `Alter`, `TimeSeries`, `Toggle` | Event/action semantics and time-series input mapping. | Untriaged | +| `Summary`, `Output` | Reporting configuration; define output-channel mapping after model compatibility. | Untriaged | + +## Numerical-regression policy + +- Keep minimized ANDES-derived input cases and captured reference results in + `test/test_files/andes_tests/`. +- The regular C++ tests must not execute ANDES. Refresh a reference only by + running the documented ANDES case deliberately and reviewing the diff. +- Each reference file stores its tolerance. Start with `1e-6` p.u. for + well-conditioned power-flow comparisons; relax only with a documented + solver-conditioning reason. +- Compare selected voltages, angles, injections, and flows. A successful + solve by itself is not sufficient evidence of compatibility. + +## Current reference cases + +| Case | Coverage | +| -------------------------- | ----------------------------------------------------------------------- | +| `andes_kundur_vsc_pflow` | 10-bus AC network, DC resistor, PQ/VQ VSC controls, and AC/DC coupling. | +| `andes_two_bus_pflow` | Minimal AC Slack/PQ/Line power flow. | +| `andes_vsc_resistor_pflow` | Minimal AC/DC `VSCShunt` plus DC resistance power flow. | diff --git a/src/fileInput/CMakeLists.txt b/src/fileInput/CMakeLists.txt index a079deaf8..3e9fbc325 100644 --- a/src/fileInput/CMakeLists.txt +++ b/src/fileInput/CMakeLists.txt @@ -42,6 +42,7 @@ set(formattedfileInput_sources gridReadPSAT.cpp readMatDyn.cpp gridDynReadCSV.cpp + gridReadAndes.cpp readMatlabData.cpp loadGDZ.cpp ) diff --git a/src/fileInput/fileInput.cpp b/src/fileInput/fileInput.cpp index 8f1d21b51..a716a21c7 100644 --- a/src/fileInput/fileInput.cpp +++ b/src/fileInput/fileInput.cpp @@ -11,6 +11,7 @@ #include "formatInterpreters/jsonReaderElement.h" #include "formatInterpreters/yamlReaderElement.h" #include "gmlc/utilities/stringOps.h" +#include "gridReadAndes.h" #include "griddyn/GridDynSimulation.h" #include "readElement.h" #include "readElementFile.h" @@ -153,13 +154,13 @@ void loadFile(std::unique_ptr& gds, } std::unique_ptr readSimXMLFile(const std::string& fileName, - ReaderInfo* ReaderInfoPtr) + ReaderInfo* readerInfoPtr) { if (!std::filesystem::exists(fileName)) { return nullptr; } return std::unique_ptr(static_cast( - loadElementFile(nullptr, fileName, ReaderInfoPtr))); + loadElementFile(nullptr, fileName, readerInfoPtr))); } void loadFile(CoreObject* parentObject, @@ -201,7 +202,9 @@ void loadFile(CoreObject* parentObject, } else if (ext == "epc") { loadEpc(parentObject, fileName, *readerInf); } else if (ext == "json") { - loadElementFile(parentObject, fileName, readerInf); + if (!loadAndesJson(parentObject, fileName)) { + loadElementFile(parentObject, fileName, readerInf); + } } else if ((ext == "yaml") || (ext == "yml")) { loadElementFile(parentObject, fileName, readerInf); } else if (ext == "gdz") { // gridDyn Zipped file diff --git a/src/fileInput/fileInput.h b/src/fileInput/fileInput.h index b1eb5b41b..21d46c01e 100644 --- a/src/fileInput/fileInput.h +++ b/src/fileInput/fileInput.h @@ -74,7 +74,7 @@ enum ReaderFlags { }; std::unique_ptr readSimXMLFile(const std::string& fileName, - ReaderInfo* ReaderInfoPtr = nullptr); + ReaderInfo* readerInfoPtr = nullptr); void addFlags(BasicReaderInfo& bri, const std::string& flags); diff --git a/src/fileInput/gridReadAndes.cpp b/src/fileInput/gridReadAndes.cpp new file mode 100644 index 000000000..f4495295b --- /dev/null +++ b/src/fileInput/gridReadAndes.cpp @@ -0,0 +1,277 @@ +/* + * Copyright (c) 2014-2026, Lawrence Livermore National Security + * See the top-level NOTICE for additional details. All rights reserved. + * SPDX-License-Identifier: BSD-3-Clause + */ + +#include "gridReadAndes.h" + +#include "core/CoreExceptions.h" +#include "griddyn/Generator.h" +#include "griddyn/GridBus.h" +#include "griddyn/Link.h" +#include "griddyn/links/AcLine.h" +#include "griddyn/links/DcLink.h" +#include "griddyn/links/VSCShunt.h" +#include "griddyn/loads/ZipLoad.h" +#include "griddyn/primary/AcBus.h" +#include "griddyn/primary/DcBus.h" +#include "nlohmann/json.hpp" +#include +#include +#include +#include +#include +#include + +namespace griddyn { +namespace { + using Json = nlohmann::json; + + std::string indexKey(const Json& record, std::string_view field = "idx") + { + if (!record.contains(field) || record[field].is_null()) { + return {}; + } + if (record[field].is_string()) { + return record[field].get(); + } + return record[field].dump(); + } + + std::string objectName(const Json& record, std::string_view fallback) + { + if (record.contains("name") && record["name"].is_string() && !record["name"].empty()) { + return record["name"].get(); + } + return std::string{fallback} + "_" + indexKey(record); + } + + double number(const Json& record, std::string_view field, double defaultValue = 0.0) + { + if (!record.contains(field) || record[field].is_null()) { + return defaultValue; + } + return record[field].get(); + } + + template + void setIfPresent(Object* object, + const Json& record, + std::string_view source, + std::string_view target) + { + if (record.contains(source) && !record[source].is_null()) { + object->set(target, number(record, source)); + } + } +} // namespace + +bool loadAndesJson(CoreObject* parentObject, const std::string& fileName) +{ + std::ifstream input(fileName); + if (!input.is_open()) { + return false; + } + + Json document; + try { + input >> document; + } + catch (const Json::parse_error&) { + return false; + } + + // Node is the distinguishing ANDES DC-topology section. Do not claim + // ordinary GridDyn JSON files that happen to contain a similarly named field. + if (!document.is_object() || !document.contains("Node") || !document["Node"].is_array()) { + return false; + } + if (parentObject == nullptr) { + throw(InvalidParameterValue("ANDES JSON import requires a simulation parent")); + } + + std::unordered_map acBuses; + std::unordered_map acBaseVoltages; + if (document.contains("Bus") && document["Bus"].is_array()) { + for (const auto& record : document["Bus"]) { + auto* bus = new AcBus(objectName(record, "Bus")); + setIfPresent(bus, record, "Vn", "basevoltage"); + setIfPresent(bus, record, "v0", "voltage"); + setIfPresent(bus, record, "a0", "angle"); + parentObject->add(bus); + acBuses.emplace(indexKey(record), bus); + acBaseVoltages.emplace(indexKey(record), number(record, "Vn", 0.0)); + } + } + + // ANDES supplies static injections separately from its AC bus records. + // Map those power-flow objects before adding the network branches. + if (document.contains("PQ") && document["PQ"].is_array()) { + for (const auto& record : document["PQ"]) { + const auto bus = acBuses.find(indexKey(record, "bus")); + if (bus == acBuses.end()) { + continue; + } + auto* load = + new ZipLoad(number(record, "p0"), number(record, "q0"), objectName(record, "PQ")); + bus->second->add(load); + } + } + if (document.contains("PV") && document["PV"].is_array()) { + for (const auto& record : document["PV"]) { + const auto bus = acBuses.find(indexKey(record, "bus")); + if (bus == acBuses.end()) { + continue; + } + bus->second->set("type", "pv"); + setIfPresent(bus->second, record, "v0", "vtarget"); + auto* generator = new Generator(objectName(record, "PV")); + generator->set("p", number(record, "p0")); + bus->second->add(generator); + } + } + if (document.contains("Slack") && document["Slack"].is_array()) { + for (const auto& record : document["Slack"]) { + const auto bus = acBuses.find(indexKey(record, "bus")); + if (bus == acBuses.end()) { + continue; + } + bus->second->set("type", "swing"); + setIfPresent(bus->second, record, "v0", "vtarget"); + setIfPresent(bus->second, record, "a0", "atarget"); + auto* generator = new Generator(objectName(record, "Slack")); + generator->set("p", number(record, "p0")); + generator->set("q", number(record, "q0")); + bus->second->add(generator); + } + } + if (document.contains("Line") && document["Line"].is_array()) { + for (const auto& record : document["Line"]) { + const auto first = acBuses.find(indexKey(record, "bus1")); + const auto second = acBuses.find(indexKey(record, "bus2")); + if ((first == acBuses.end()) || (second == acBuses.end())) { + continue; + } + auto* line = new AcLine(objectName(record, "Line")); + setIfPresent(line, record, "r", "r"); + setIfPresent(line, record, "x", "x"); + setIfPresent(line, record, "b", "b"); + setIfPresent(line, record, "tap", "tap"); + setIfPresent(line, record, "phi", "tapangle"); + line->updateBus(first->second, 1); + line->updateBus(second->second, 2); + if (number(record, "u", 1.0) == 0.0) { + line->disable(); + } + parentObject->add(line); + } + } + + std::unordered_map dcBuses; + for (const auto& record : document["Node"]) { + auto* bus = new DcBus(objectName(record, "Node")); + setIfPresent(bus, record, "Vdcn", "basevoltage"); + setIfPresent(bus, record, "v0", "voltage"); + parentObject->add(bus); + dcBuses.emplace(indexKey(record), bus); + } + + if (document.contains("Ground") && document["Ground"].is_array()) { + for (const auto& record : document["Ground"]) { + const auto node = indexKey(record, "node"); + const auto found = dcBuses.find(node); + if (found != dcBuses.end()) { + found->second->set("type", "swing"); + setIfPresent(found->second, record, "voltage", "voltage"); + } + } + } + + const std::array, 8> branchModels{{ + {"R", "r"}, + {"L", "l"}, + {"C", "c"}, + {"RLs", "rls"}, + {"RCp", "rcp"}, + {"RLCp", "rlcp"}, + {"RCs", "rcs"}, + {"RLCs", "rlcs"}, + }}; + for (const auto& [section, model] : branchModels) { + if (!document.contains(section) || !document[section].is_array()) { + continue; + } + for (const auto& record : document[section]) { + const auto first = dcBuses.find(indexKey(record, "node1")); + const auto second = dcBuses.find(indexKey(record, "node2")); + if ((first == dcBuses.end()) || (second == dcBuses.end())) { + continue; + } + auto* branch = new links::DcLink(objectName(record, section)); + branch->set("model", model); + branch->set("andes_current_balance", 1.0); + setIfPresent(branch, record, "R", "r"); + setIfPresent(branch, record, "L", "l"); + setIfPresent(branch, record, "C", "c"); + // ANDES' R model has no implicit inductance; GridDyn's historical + // DcLink default does, so explicitly remove it during conversion. + if (model == "r") { + branch->set("l", 0.0); + } + branch->updateBus(first->second, 1); + branch->updateBus(second->second, 2); + parentObject->add(branch); + } + } + + if (document.contains("VSCShunt") && document["VSCShunt"].is_array()) { + for (const auto& record : document["VSCShunt"]) { + const auto acBus = acBuses.find(indexKey(record, "bus")); + const auto dcBus = dcBuses.find(indexKey(record, "node1")); + const auto dcReference = dcBuses.find(indexKey(record, "node2")); + if ((acBus == acBuses.end()) || (dcBus == dcBuses.end()) || + (dcReference == dcBuses.end())) { + continue; + } + auto* converter = new links::VSCShunt(objectName(record, "VSCShunt")); + converter->set("andes_current_balance", 1.0); + // In Andes, rsh and xsh are impedance parameters (z=True). They + // are converted to the AC bus base using (VSC Vn / Bus Vn)^2. + // GridDyn stores the values directly on the bus base, so carry + // out that conversion while importing the original JSON values. + const auto baseVoltage = acBaseVoltages.find(indexKey(record, "bus")); + const auto busVn = (baseVoltage != acBaseVoltages.end()) ? baseVoltage->second : 0.0; + const auto vscVn = number(record, "Vn", busVn); + const auto impedanceScale = + (busVn > 0.0 && vscVn > 0.0) ? std::pow(vscVn / busVn, 2) : 1.0; + if (record.contains("rsh") && !record["rsh"].is_null()) { + converter->set("r", number(record, "rsh") * impedanceScale); + } + if (record.contains("xsh") && !record["xsh"].is_null()) { + converter->set("x", number(record, "xsh") * impedanceScale); + } + setIfPresent(converter, record, "control", "control"); + setIfPresent(converter, record, "v0", "v0"); + setIfPresent(converter, record, "p0", "p0"); + setIfPresent(converter, record, "q0", "q0"); + setIfPresent(converter, record, "vdc0", "vdc0"); + setIfPresent(converter, record, "k0", "k0"); + setIfPresent(converter, record, "k1", "k1"); + setIfPresent(converter, record, "k2", "k2"); + setIfPresent(converter, record, "droop", "droop"); + setIfPresent(converter, record, "K", "k"); + setIfPresent(converter, record, "vhigh", "vhigh"); + setIfPresent(converter, record, "vlow", "vlow"); + setIfPresent(converter, record, "vshmax", "vshmax"); + setIfPresent(converter, record, "vshmin", "vshmin"); + setIfPresent(converter, record, "Ishmax", "ishmax"); + converter->updateBus(acBus->second, 1); + converter->updateBus(dcBus->second, 2); + converter->updateBus(dcReference->second, 3); + parentObject->add(converter); + } + } + return true; +} +} // namespace griddyn diff --git a/src/fileInput/gridReadAndes.h b/src/fileInput/gridReadAndes.h new file mode 100644 index 000000000..64a09c6e5 --- /dev/null +++ b/src/fileInput/gridReadAndes.h @@ -0,0 +1,20 @@ +/* + * Copyright (c) 2014-2026, Lawrence Livermore National Security + * See the top-level NOTICE for additional details. All rights reserved. + * SPDX-License-Identifier: BSD-3-Clause + */ + +#pragma once + +#include + +namespace griddyn { +class CoreObject; + +/** Load the DC-model sections used by an ANDES JSON case. + * + * Returns false when the file is not an ANDES case, allowing callers to use + * the normal GridDyn JSON reader for native files. + */ +bool loadAndesJson(CoreObject* parentObject, const std::string& fileName); +} // namespace griddyn diff --git a/src/griddyn/CMakeLists.txt b/src/griddyn/CMakeLists.txt index 1996f368b..2744ce563 100644 --- a/src/griddyn/CMakeLists.txt +++ b/src/griddyn/CMakeLists.txt @@ -202,6 +202,7 @@ set(load_sources set(link_headers links/DcLink.h links/AcDcConverter.h + links/VSCShunt.h links/AdjustableTransformer.h Link.h links/Subsystem.h @@ -218,6 +219,7 @@ set(link_sources links/AdjustableTransformer.cpp links/DcLink.cpp links/AcDcConverter.cpp + links/VSCShunt.cpp links/Subsystem.cpp links/Hvdc.cpp links/ZBreaker.cpp diff --git a/src/griddyn/links/DcLink.cpp b/src/griddyn/links/DcLink.cpp index 22797303e..f33a544d3 100644 --- a/src/griddyn/links/DcLink.cpp +++ b/src/griddyn/links/DcLink.cpp @@ -11,8 +11,10 @@ #include "../primary/DcBus.h" #include "core/CoreExceptions.h" #include "core/CoreObjectTemplates.hpp" +#include "gmlc/utilities/stringOps.h" #include "gmlc/utilities/vectorOps.hpp" #include +#include #include #include @@ -41,10 +43,32 @@ CoreObject* DcLink::clone(CoreObject* obj) const } nobj->r = r; nobj->x = x; + nobj->c = c; + nobj->topology = topology; + nobj->currentBalance = currentBalance; return nobj; } +bool DcLink::hasCapacitor() const +{ + return topology != Topology::RL; +} + +bool DcLink::isOpenInPowerFlow() const +{ + return (topology == Topology::C) || (topology == Topology::RC_SERIES) || + (topology == Topology::RLC_SERIES); +} + +count_t DcLink::currentStateOffset(const SolverMode& sMode) const +{ + if (hasCapacitor()) { + return offsets.getAlgOffset(sMode); + } + return isDynamic(sMode) ? offsets.getDiffOffset(sMode) : offsets.getAlgOffset(sMode); +} + void DcLink::timestep(CoreTime /*time*/, const IOdata& /*inputs*/, const SolverMode& /*sMode*/) { if (!isEnabled()) { @@ -87,7 +111,30 @@ double DcLink::getMaxTransfer() const // set properties void DcLink::set(std::string_view param, std::string_view val) { - Link::set(param, val); + if ((param == "model") || (param == "topology")) { + const auto model = gmlc::utilities::convertToLowerCase(val); + if ((model == "r") || (model == "l") || (model == "rls") || (model == "rl")) { + topology = Topology::RL; + } else if (model == "c") { + topology = Topology::C; + } else if (model == "rcp") { + topology = Topology::RC_PARALLEL; + } else if (model == "rlcp") { + topology = Topology::RLC_PARALLEL; + } else if (model == "rcs") { + topology = Topology::RC_SERIES; + } else if (model == "rlcs") { + topology = Topology::RLC_SERIES; + } else { + throw(InvalidParameterValue(val)); + } + } else if ((param == "andes_current_balance") || (param == "current_balance")) { + const auto mode = gmlc::utilities::convertToLowerCase(val); + currentBalance = + (mode == "true") || (mode == "1") || (mode == "current") || (mode == "andes"); + } else { + Link::set(param, val); + } } void DcLink::set(std::string_view param, double val, unit unitType) { @@ -96,6 +143,10 @@ void DcLink::set(std::string_view param, double val, unit unitType) } else if ((param == "l") || (param == "x") || (param == "ldc")) { x = val; // set line admittance + } else if ((param == "c") || (param == "cdc")) { + c = val; + } else if ((param == "andes_current_balance") || (param == "current_balance")) { + currentBalance = (val > 0.5); } else { Link::set(param, val, unitType); } @@ -123,13 +174,21 @@ void DcLink::pFlowObjectInitializeB() void DcLink::dynObjectInitializeA(CoreTime /*time0*/, std::uint32_t /*flags*/) { - m_dstate_dt.resize(1); - m_state.resize(1); + std::uint32_t differentialCount = 0U; + if (topology == Topology::RLC_PARALLEL || topology == Topology::RLC_SERIES) { + differentialCount = 2U; + } else if (hasCapacitor() || (x != 0.0)) { + differentialCount = 1U; + } + m_dstate_dt.assign(differentialCount, 0.0); + m_state.assign(differentialCount, 0.0); updateLocalCache(); - if (x != 0) { + if (!hasCapacitor() && (x != 0)) { m_state[0] = Idc; - m_dstate_dt[0] = 0.0; - } else { + } else if (hasCapacitor()) { + // The capacitor voltage is always the first differential state. + m_state[(topology == Topology::RLC_PARALLEL || topology == Topology::RLC_SERIES) ? 1 : 0] = + linkInfo.v1 - linkInfo.v2; } } @@ -137,12 +196,17 @@ StateSizes DcLink::localStateSizes(const SolverMode& sMode) const { StateSizes localSS; if (isDynamic(sMode)) { - if (x != 0.0) { + if (hasCapacitor()) { if (!isAlgebraicOnly(sMode)) { - localSS.diffSize = 1; + localSS.diffSize = + (topology == Topology::RLC_PARALLEL || topology == Topology::RLC_SERIES) ? 2 : + 1; } + localSS.algSize = 1; + } else if ((x != 0.0) && !isAlgebraicOnly(sMode)) { + localSS.diffSize = 1; } - } else if (r <= 0.0) { // superconducting + } else if (!isOpenInPowerFlow() && (r <= 0.0)) { // superconducting localSS.algSize = 1; } return localSS; @@ -152,12 +216,14 @@ count_t DcLink::localJacobianCount(const SolverMode& sMode) const { count_t jacCount = 0; if (isDynamic(sMode)) { - if (x != 0.0) { - if (!isAlgebraicOnly(sMode)) { - jacCount = 3; - } + if (hasCapacitor()) { + const auto diff = + (topology == Topology::RLC_PARALLEL || topology == Topology::RLC_SERIES) ? 2 : 1; + jacCount = static_cast(4 + (2 * diff)); + } else if ((x != 0.0) && !isAlgebraicOnly(sMode)) { + jacCount = 3; } - } else if (r <= 0.0) { // superconducting + } else if (!isOpenInPowerFlow() && (r <= 0.0)) { // superconducting jacCount = 2; } return jacCount; @@ -175,6 +241,13 @@ void DcLink::ioPartialDerivatives(id_type_t busId, return; } + if (currentBalance) { + if (stateSize(sMode) > 0) { + return; + } + jacobian.assignCheckCol(POUT_LOCATION, inputLocs[VOLTAGE_IN_LOCATION], 1.0 / r); + return; + } if ((busId == 2) || (busId == B2->getID())) { jacobian.assignCheckCol(POUT_LOCATION, inputLocs[VOLTAGE_IN_LOCATION], -Idc); } else { @@ -195,7 +268,7 @@ void DcLink::outputPartialDerivatives(id_type_t busId, double p1v2 = 0.0; double p2v1 = 0.0; if (!isDynamic(sMode)) { - if (r > 0.0) { + if ((r > 0.0) && !isOpenInPowerFlow()) { p1v2 = -linkInfo.v1 / r; p2v1 = -linkInfo.v1 / r; } else { @@ -206,20 +279,34 @@ void DcLink::outputPartialDerivatives(id_type_t busId, // md.assign(B1Voffset, B2Voffset, Q1V2); // md.assign(B2Voffset, B1Voffset, Q2V1); + if (currentBalance) { + if ((busId == 2) || (busId == B2->getID())) { + if (stateSize(sMode) > 0) { + jacobian.assign(POUT_LOCATION, currentStateOffset(sMode), -1.0); + } else if (r > 0.0) { + const int bus1VoltageOffset = B1->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); + jacobian.assignCheckCol(POUT_LOCATION, bus1VoltageOffset, -1.0 / r); + } + } else { + if (stateSize(sMode) > 0) { + jacobian.assign(POUT_LOCATION, currentStateOffset(sMode), 1.0); + } else if (r > 0.0) { + const int bus2VoltageOffset = B2->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); + jacobian.assignCheckCol(POUT_LOCATION, bus2VoltageOffset, -1.0 / r); + } + } + return; + } if ((busId == 2) || (busId == B2->getID())) { if (stateSize(sMode) > 0) { - auto offset = - isDynamic(sMode) ? offsets.getDiffOffset(sMode) : offsets.getAlgOffset(sMode); - jacobian.assign(POUT_LOCATION, offset, -linkInfo.v2); + jacobian.assign(POUT_LOCATION, currentStateOffset(sMode), -linkInfo.v2); } else { const int bus1VoltageOffset = B1->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); jacobian.assignCheckCol(POUT_LOCATION, bus1VoltageOffset, p2v1); } } else { if (stateSize(sMode) > 0) { - auto offset = - isDynamic(sMode) ? offsets.getDiffOffset(sMode) : offsets.getAlgOffset(sMode); - jacobian.assign(POUT_LOCATION, offset, linkInfo.v1); + jacobian.assign(POUT_LOCATION, currentStateOffset(sMode), linkInfo.v1); } else { const int bus2VoltageOffset = B2->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); jacobian.assignCheckCol(POUT_LOCATION, bus2VoltageOffset, p1v2); @@ -238,7 +325,46 @@ void DcLink::jacobianElements(const IOdata& /*inputs*/, const IOlocs& /*inputLocs*/, const SolverMode& sMode) { - if (stateSize(sMode) > 0) { + if (hasCapacitor() && isDynamic(sMode)) { + const int bus1VoltageOffset = B1->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); + const int bus2VoltageOffset = B2->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); + const auto currentOffset = offsets.getAlgOffset(sMode); + const auto firstDiffOffset = offsets.getDiffOffset(sMode); + const auto capacitorOffset = firstDiffOffset + + ((topology == Topology::RLC_PARALLEL || topology == Topology::RLC_SERIES) ? 1 : 0); + + // Capacitor voltage constraint: vC = v1 - v2, with Idc as a branch current. + jacobian.assign(currentOffset, capacitorOffset, 1.0); + jacobian.assignCheckCol(currentOffset, bus1VoltageOffset, -1.0); + jacobian.assignCheckCol(currentOffset, bus2VoltageOffset, 1.0); + + if (topology == Topology::C) { + jacobian.assign(capacitorOffset, capacitorOffset, -stateData.cj); + jacobian.assign(capacitorOffset, currentOffset, -1.0 / c); + } else if (topology == Topology::RC_PARALLEL) { + jacobian.assign(capacitorOffset, capacitorOffset, (1.0 / (r * c)) - stateData.cj); + jacobian.assign(capacitorOffset, currentOffset, -1.0 / c); + } else if (topology == Topology::RLC_PARALLEL) { + jacobian.assign(firstDiffOffset, capacitorOffset, 1.0 / x); + jacobian.assign(firstDiffOffset, firstDiffOffset, -stateData.cj); + jacobian.assign(capacitorOffset, firstDiffOffset, 1.0 / c); + jacobian.assign(capacitorOffset, capacitorOffset, (1.0 / (r * c)) - stateData.cj); + jacobian.assign(capacitorOffset, currentOffset, -1.0 / c); + } else if (topology == Topology::RC_SERIES) { + jacobian.assign(currentOffset, currentOffset, -r); + jacobian.assign(capacitorOffset, capacitorOffset, -stateData.cj); + jacobian.assign(capacitorOffset, currentOffset, -1.0 / c); + } else { // RLC series + jacobian.assign(currentOffset, firstDiffOffset, -1.0); + jacobian.assign(currentOffset, currentOffset, -1.0); + jacobian.assign(firstDiffOffset, firstDiffOffset, (-r / x) - stateData.cj); + jacobian.assign(firstDiffOffset, capacitorOffset, -1.0 / x); + jacobian.assignCheckCol(firstDiffOffset, bus1VoltageOffset, 1.0 / x); + jacobian.assignCheckCol(firstDiffOffset, bus2VoltageOffset, -1.0 / x); + jacobian.assign(capacitorOffset, firstDiffOffset, 1.0 / c); + jacobian.assign(capacitorOffset, capacitorOffset, -stateData.cj); + } + } else if (stateSize(sMode) > 0) { const int bus1VoltageOffset = B1->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); const int bus2VoltageOffset = B2->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); updateLocalCache(noInputs, stateData, sMode); @@ -266,7 +392,38 @@ void DcLink::residual(const IOdata& inputs, double resid[], const SolverMode& sMode) { - if (stateSize(sMode) > 0) { + if (hasCapacitor() && isDynamic(sMode)) { + updateLocalCache(inputs, stateData, sMode); + const auto currentOffset = offsets.getAlgOffset(sMode); + const auto firstDiffOffset = offsets.getDiffOffset(sMode); + const auto capacitorOffset = firstDiffOffset + + ((topology == Topology::RLC_PARALLEL || topology == Topology::RLC_SERIES) ? 1 : 0); + const auto capacitorVoltage = stateData.state[capacitorOffset]; + const auto current = stateData.state[currentOffset]; + resid[currentOffset] = capacitorVoltage - (linkInfo.v1 - linkInfo.v2); + + if (topology == Topology::C) { + resid[capacitorOffset] = (-current / c) - stateData.dstate_dt[capacitorOffset]; + } else if (topology == Topology::RC_PARALLEL) { + resid[capacitorOffset] = + (-(current - (capacitorVoltage / r)) / c) - stateData.dstate_dt[capacitorOffset]; + } else if (topology == Topology::RLC_PARALLEL) { + const auto inductorCurrent = stateData.state[firstDiffOffset]; + resid[firstDiffOffset] = (capacitorVoltage / x) - stateData.dstate_dt[firstDiffOffset]; + resid[capacitorOffset] = (-(current - (capacitorVoltage / r) - inductorCurrent) / c) - + stateData.dstate_dt[capacitorOffset]; + } else if (topology == Topology::RC_SERIES) { + resid[currentOffset] -= r * current; + resid[capacitorOffset] = (-current / c) - stateData.dstate_dt[capacitorOffset]; + } else { + const auto inductorCurrent = stateData.state[firstDiffOffset]; + resid[currentOffset] = -inductorCurrent - current; + resid[firstDiffOffset] = + ((linkInfo.v1 - linkInfo.v2 - (r * inductorCurrent) - capacitorVoltage) / x) - + stateData.dstate_dt[firstDiffOffset]; + resid[capacitorOffset] = (inductorCurrent / c) - stateData.dstate_dt[capacitorOffset]; + } + } else if (stateSize(sMode) > 0) { updateLocalCache(inputs, stateData, sMode); if (isDynamic(sMode)) { auto offset = offsets.getDiffOffset(sMode); @@ -289,7 +446,16 @@ void DcLink::setState(CoreTime time, const double dstateDt[], const SolverMode& sMode) { - if (stateSize(sMode) > 0) { + if (hasCapacitor() && isDynamic(sMode)) { + const auto firstDiffOffset = offsets.getDiffOffset(sMode); + const auto differentialCount = + (topology == Topology::RLC_PARALLEL || topology == Topology::RLC_SERIES) ? 2U : 1U; + for (unsigned int ii = 0; ii < differentialCount; ++ii) { + m_state[ii] = state[firstDiffOffset + ii]; + m_dstate_dt[ii] = dstateDt[firstDiffOffset + ii]; + } + Idc = state[offsets.getAlgOffset(sMode)]; + } else if (stateSize(sMode) > 0) { if (isDynamic(sMode)) { auto offset = offsets.getDiffOffset(sMode); m_state[0] = state[offset]; @@ -308,7 +474,14 @@ void DcLink::guessState(const CoreTime /*time*/, double dstateDt[], const SolverMode& sMode) { - if (stateSize(sMode) > 0) { + if (hasCapacitor() && isDynamic(sMode)) { + const auto firstDiffOffset = offsets.getDiffOffset(sMode); + for (unsigned int ii = 0; ii < m_state.size(); ++ii) { + state[firstDiffOffset + ii] = m_state[ii]; + dstateDt[firstDiffOffset + ii] = m_dstate_dt[ii]; + } + state[offsets.getAlgOffset(sMode)] = Idc; + } else if (stateSize(sMode) > 0) { if (isDynamic(sMode)) { auto offset = offsets.getDiffOffset(sMode); state[offset] = m_state[0]; @@ -324,7 +497,17 @@ void DcLink::getStateName(stringVec& stNames, const SolverMode& sMode, const std::string& prefix) const { - if (stateSize(sMode) > 0) { + if (hasCapacitor() && isDynamic(sMode)) { + const std::string prefix2 = prefix + getName() + ':'; + const auto firstDiffOffset = offsets.getDiffOffset(sMode); + if (topology == Topology::RLC_PARALLEL || topology == Topology::RLC_SERIES) { + stNames[firstDiffOffset] = prefix2 + "il"; + stNames[firstDiffOffset + 1] = prefix2 + "vc"; + } else { + stNames[firstDiffOffset] = prefix2 + "vc"; + } + stNames[offsets.getAlgOffset(sMode)] = prefix2 + "idc"; + } else if (stateSize(sMode) > 0) { const std::string prefix2 = prefix + getName() + ':'; auto offset = (isDynamic(sMode)) ? offsets.getDiffOffset(sMode) : offsets.getAlgOffset(sMode); @@ -348,11 +531,11 @@ void DcLink::updateLocalCache(const IOdata& /*inputs*/, linkInfo.v1 = B1->getVoltage(stateData.state, sMode); linkInfo.v2 = B2->getVoltage(stateData.state, sMode); if (stateSize(sMode) > 0) { - auto offset = - (isDynamic(sMode)) ? offsets.getDiffOffset(sMode) : offsets.getAlgOffset(sMode); - Idc = stateData.state[offset]; + Idc = stateData.state[currentStateOffset(sMode)]; } else { - if (r > 0) { + if (isOpenInPowerFlow()) { + Idc = 0.0; + } else if (r > 0) { Idc = (linkInfo.v1 - linkInfo.v2) / r; // Q2 = P2*sqrt(k3sq2*k3sq2 - gamma*gamma); @@ -360,8 +543,8 @@ void DcLink::updateLocalCache(const IOdata& /*inputs*/, Idc = Pset / linkInfo.v1; } } - linkFlows.P1 = linkInfo.v1 * Idc; - linkFlows.P2 = -linkInfo.v2 * Idc; + linkFlows.P1 = currentBalance ? Idc : linkInfo.v1 * Idc; + linkFlows.P2 = currentBalance ? -Idc : -linkInfo.v2 * Idc; } void DcLink::updateLocalCache() @@ -371,8 +554,8 @@ void DcLink::updateLocalCache() if (isEnabled()) { linkInfo.v1 = B1->getVoltage(); linkInfo.v2 = B2->getVoltage(); - linkFlows.P1 = linkInfo.v1 * Idc; - linkFlows.P2 = -linkInfo.v2 * Idc; + linkFlows.P1 = currentBalance ? Idc : linkInfo.v1 * Idc; + linkFlows.P2 = currentBalance ? -Idc : -linkInfo.v2 * Idc; } } diff --git a/src/griddyn/links/DcLink.h b/src/griddyn/links/DcLink.h index 6fcd1ae91..50f5eb8ac 100644 --- a/src/griddyn/links/DcLink.h +++ b/src/griddyn/links/DcLink.h @@ -14,6 +14,8 @@ namespace griddyn::links { */ class DcLink: public Link { public: + /** DC branch topologies compatible with the ANDES dc.rlc models. */ + enum class Topology { RL, C, RC_PARALLEL, RLC_PARALLEL, RC_SERIES, RLC_SERIES }; /* enum DcLinkFlags { FIXED_TARGET_POWER = OBJECT_FLAG5, @@ -22,6 +24,14 @@ class DcLink: public Link { double Idc = 0; //!< [puA] storage for DC current double r = 0; //!< [puOhm] the dc resistance double x = 0.0001; //!< [puOhm] the dc inductance + double c = 0.0; //!< [puF] the dc capacitance + Topology topology = Topology::RL; + bool currentBalance = false; //!< use ANDES' nodal-current convention + + bool hasCapacitor() const; + bool isOpenInPowerFlow() const; + count_t currentStateOffset(const SolverMode& sMode) const; + public: DcLink(const std::string& objName = "dclink_$"); DcLink(double resistancePu, double reactancePu, const std::string& objName = "dclink_$"); diff --git a/src/griddyn/links/Link.cpp b/src/griddyn/links/Link.cpp index 8dd2773c2..5a2862794 100644 --- a/src/griddyn/links/Link.cpp +++ b/src/griddyn/links/Link.cpp @@ -12,6 +12,7 @@ #include "AcLine.h" #include "AdjustableTransformer.h" #include "DcLink.h" +#include "VSCShunt.h" #include "core/CoreExceptions.h" #include "core/CoreObjectTemplates.hpp" #include "core/ObjectFactoryTemplates.hpp" @@ -62,6 +63,8 @@ namespace links { AcDcConverter::Mode::INVERTER); static ChildTypeFactory gAcdc("link", std::to_array({"acdc", "acdcconverter", "dcconverter"})); + static ChildTypeFactory + gVscShunt("link", std::to_array({"vscshunt", "vsc_shunt"})); } // namespace links std::atomic Link::linkCount(0); // helper defines to have things make more sense diff --git a/src/griddyn/links/VSCShunt.cpp b/src/griddyn/links/VSCShunt.cpp new file mode 100644 index 000000000..015dc2d85 --- /dev/null +++ b/src/griddyn/links/VSCShunt.cpp @@ -0,0 +1,496 @@ +/* + * Copyright (c) 2014-2026, Lawrence Livermore National Security + * See the top-level NOTICE for additional details. All rights reserved. + * SPDX-License-Identifier: BSD-3-Clause + */ + +#include "VSCShunt.h" + +#include "../GridBus.h" +#include "../primary/DcBus.h" +#include "core/CoreExceptions.h" +#include "core/CoreObjectTemplates.hpp" +#include "gmlc/utilities/stringOps.h" +#include +#include +#include + +namespace griddyn::links { +using units::unit; + +namespace { + constexpr index_t ashIndex = 0; + constexpr index_t vshIndex = 1; + constexpr index_t pshIndex = 2; + constexpr index_t qshIndex = 3; + constexpr index_t pdcIndex = 4; +} // namespace + +VSCShunt::VSCShunt(const std::string& objName): AcDcConverter(objName) {} + +CoreObject* VSCShunt::clone(CoreObject* obj) const +{ + auto* nobj = cloneBase(this, obj); + if (nobj == nullptr) { + return obj; + } + nobj->control = control; + nobj->v0 = v0; + nobj->p0 = p0; + nobj->q0 = q0; + nobj->vdc0 = vdc0; + nobj->k0 = k0; + nobj->k1 = k1; + nobj->k2 = k2; + nobj->droop = droop; + nobj->droopK = droopK; + nobj->vhigh = vhigh; + nobj->vlow = vlow; + nobj->vshmax = vshmax; + nobj->vshmin = vshmin; + nobj->ishmax = ishmax; + nobj->currentBalance = currentBalance; + nobj->ash = ash; + nobj->vsh = vsh; + nobj->psh = psh; + nobj->qsh = qsh; + nobj->pdc = pdc; + return nobj; +} + +bool VSCShunt::isConnected() const +{ + return AcDcConverter::isConnected() && (dcReference != nullptr); +} + +GridBus* VSCShunt::getBus(index_t busInd) const +{ + if (busInd == 3) { + return dcReference; + } + return AcDcConverter::getBus(busInd); +} + +void VSCShunt::followNetwork(int network, std::queue& stk) +{ + // AcDcConverter intentionally clears NETWORK_CONNECTED so conventional + // two-terminal converters do not merge AC and DC islands. VSCShunt has + // an explicit DC reference terminal and algebraic equations spanning all + // three terminals, so its complete electrical island must be traversed. + if (!isConnected()) { + return; + } + for (auto* bus : {B1, B2, static_cast(dcReference)}) { + if ((bus != nullptr) && (bus->Network != network)) { + stk.push(bus); + } + } +} + +void VSCShunt::updateBus(GridBus* bus, index_t busnumber) +{ + if (busnumber == 3) { + auto* dcBus = dynamic_cast(bus); + if (dcBus == nullptr) { + throw(UnrecognizedObjectException(this)); + } + if (dcReference != nullptr) { + dcReference->remove(this); + } + dcReference = dcBus; + dcReference->add(this); + return; + } + AcDcConverter::updateBus(bus, busnumber); +} + +void VSCShunt::set(std::string_view param, std::string_view val) +{ + const auto lower = gmlc::utilities::convertToLowerCase(val); + if (param == "control") { + if (lower == "pq") { + control = Control::PQ; + } else if (lower == "pv") { + control = Control::PV; + } else if ((lower == "vq") || (lower == "v_q")) { + control = Control::VQ; + } else if ((lower == "vv") || (lower == "v_v")) { + control = Control::VV; + } else { + throw(InvalidParameterValue(val)); + } + return; + } + AcDcConverter::set(param, val); +} + +void VSCShunt::set(std::string_view param, double val, unit unitType) +{ + if (param == "control") { + const auto controlValue = static_cast(val); + if ((controlValue < 0) || (controlValue > 3)) { + throw(InvalidParameterValue("control")); + } + control = static_cast(controlValue); + } else if (param == "v0") { + v0 = val; + } else if (param == "p0") { + p0 = val; + } else if (param == "q0") { + q0 = val; + } else if (param == "vdc0") { + vdc0 = val; + } else if (param == "k0") { + k0 = val; + } else if (param == "k1") { + k1 = val; + } else if (param == "k2") { + k2 = val; + } else if (param == "droop") { + droop = val; + } else if (param == "k") { + droopK = val; + } else if (param == "vhigh") { + vhigh = val; + } else if (param == "vlow") { + vlow = val; + } else if (param == "vshmax") { + vshmax = val; + } else if (param == "vshmin") { + vshmin = val; + } else if ((param == "ishmax") || (param == "imax")) { + ishmax = val; + } else if ((param == "andes_current_balance") || (param == "current_balance")) { + currentBalance = (val > 0.5); + } else { + AcDcConverter::set(param, val, unitType); + } +} + +void VSCShunt::pFlowObjectInitializeA(CoreTime time0, std::uint32_t flags) +{ + Link::pFlowObjectInitializeA(time0, flags); + ash = B1->getAngle(); + vsh = v0; + psh = p0; + qsh = q0; + pdc = 0.0; + offsets.local().local.algSize = 5; + offsets.local().local.jacSize = localJacobianCount(cPflowSolverMode); + updateLocalCache(); +} + +void VSCShunt::dynObjectInitializeA(CoreTime time0, std::uint32_t flags) +{ + // ANDES' VSCShunt is a static power-flow component. Retain its algebraic + // equations in DAE mode, without inheriting AcDcConverter's unrelated + // firing-angle controller states. + Link::dynObjectInitializeA(time0, flags); + offsets.local().local.algSize = 5; + offsets.local().local.jacSize = localJacobianCount(cDaeSolverMode); + updateLocalCache(); +} + +StateSizes VSCShunt::localStateSizes(const SolverMode& sMode) const +{ + StateSizes sizes; + if (hasAlgebraic(sMode)) { + sizes.algSize = 5; + } + return sizes; +} + +count_t VSCShunt::localJacobianCount(const SolverMode& sMode) const +{ + return hasAlgebraic(sMode) ? 19 : 0; +} + +double VSCShunt::dcVoltageDifference() const +{ + return B2->getVoltage() - static_cast(dcReference)->getVoltage(); +} + +void VSCShunt::updateFlows(double dcVoltage, double dcReferenceVoltage) +{ + linkFlows.P1 = -psh; + linkFlows.Q1 = -qsh; + const double difference = dcVoltage - dcReferenceVoltage; + if (std::abs(difference) < 1e-12) { + linkFlows.P2 = 0.0; + linkFlows.Q2 = 0.0; + return; + } + linkFlows.P2 = currentBalance ? -pdc / difference : -pdc * dcVoltage / difference; + linkFlows.Q2 = 0.0; +} + +void VSCShunt::updateLocalCache() +{ + if (!isEnabled() || !isConnected()) { + return; + } + linkInfo = {}; + linkInfo.v1 = B1->getVoltage(); + linkInfo.v2 = B2->getVoltage(); + linkInfo.theta1 = B1->getAngle(); + updateFlows(linkInfo.v2, static_cast(dcReference)->getVoltage()); +} + +void VSCShunt::updateLocalCache(const IOdata& /*inputs*/, + const StateData& stateDataValue, + const SolverMode& sMode) +{ + if (!isEnabled() || !isConnected() || !stateDataValue.updateRequired(linkInfo.seqID)) { + return; + } + linkInfo = {}; + linkInfo.seqID = stateDataValue.seqID; + linkInfo.v1 = B1->getVoltage(stateDataValue, sMode); + linkInfo.v2 = B2->getVoltage(stateDataValue, sMode); + linkInfo.theta1 = B1->getAngle(stateDataValue, sMode); + const auto offset = offsets.getAlgOffset(sMode); + ash = stateDataValue.state[offset + ashIndex]; + vsh = stateDataValue.state[offset + vshIndex]; + psh = stateDataValue.state[offset + pshIndex]; + qsh = stateDataValue.state[offset + qshIndex]; + pdc = stateDataValue.state[offset + pdcIndex]; + updateFlows(linkInfo.v2, dcReference->getVoltage(stateDataValue, sMode)); +} + +void VSCShunt::ioPartialDerivatives(id_type_t busId, + const StateData& stateDataValue, + MatrixData& matrixDataValue, + const IOlocs& inputLocs, + const SolverMode& sMode) +{ + if (!isEnabled() || (inputLocs[VOLTAGE_IN_LOCATION] == kNullLocation)) { + return; + } + updateLocalCache(noInputs, stateDataValue, sMode); + const auto difference = + B2->getVoltage(stateDataValue, sMode) - dcReference->getVoltage(stateDataValue, sMode); + if (std::abs(difference) < 1e-12) { + return; + } + if (busId == B2->getID()) { + matrixDataValue.assign(POUT_LOCATION, + inputLocs[VOLTAGE_IN_LOCATION], + currentBalance ? + pdc / (difference * difference) : + pdc * dcReference->getVoltage(stateDataValue, sMode) / + (difference * difference)); + } else if (busId == dcReference->getID()) { + matrixDataValue.assign(POUT_LOCATION, + inputLocs[VOLTAGE_IN_LOCATION], + currentBalance ? pdc / (difference * difference) : + pdc * B2->getVoltage(stateDataValue, sMode) / + (difference * difference)); + } +} + +void VSCShunt::outputPartialDerivatives(id_type_t busId, + const StateData& stateDataValue, + MatrixData& matrixDataValue, + const SolverMode& sMode) +{ + if (!isEnabled()) { + return; + } + updateLocalCache(noInputs, stateDataValue, sMode); + const auto offset = offsets.getAlgOffset(sMode); + if (busId == B1->getID()) { + matrixDataValue.assign(POUT_LOCATION, offset + pshIndex, -1.0); + matrixDataValue.assign(QOUT_LOCATION, offset + qshIndex, -1.0); + return; + } + const auto difference = + B2->getVoltage(stateDataValue, sMode) - dcReference->getVoltage(stateDataValue, sMode); + if (std::abs(difference) < 1e-12) { + return; + } + if (busId == B2->getID()) { + matrixDataValue.assign(POUT_LOCATION, + offset + pdcIndex, + currentBalance ? + -1.0 / difference : + -B2->getVoltage(stateDataValue, sMode) / difference); + } else if (busId == dcReference->getID()) { + matrixDataValue.assign(POUT_LOCATION, + offset + pdcIndex, + currentBalance ? + 1.0 / difference : + dcReference->getVoltage(stateDataValue, sMode) / difference); + } +} + +count_t VSCShunt::outputDependencyCount(index_t /*num*/, const SolverMode& sMode) const +{ + return hasAlgebraic(sMode) ? 2 : 0; +} + +void VSCShunt::jacobianElements(const IOdata& /*inputs*/, + const StateData& stateDataValue, + MatrixData& matrixDataValue, + const IOlocs& /*inputLocs*/, + const SolverMode& sMode) +{ + updateLocalCache(noInputs, stateDataValue, sMode); + const auto offset = offsets.getAlgOffset(sMode); + const auto busVoltage = B1->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); + const auto busAngle = B1->getOutputLoc(sMode, ANGLE_IN_LOCATION); + const auto dcVoltage = B2->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); + const auto dcReferenceVoltage = dcReference->getOutputLoc(sMode, VOLTAGE_IN_LOCATION); + const auto voltage = linkInfo.v1; + const auto angleDifference = linkInfo.theta1 - ash; + const auto cosine = std::cos(angleDifference); + const auto sine = std::sin(angleDifference); + const auto denominator = (r * r) + (x * x); + const auto conductance = r / denominator; + const auto susceptance = -x / denominator; + + matrixDataValue.assignCheckCol(offset + ashIndex, + busVoltage, + (2.0 * conductance * voltage) - (conductance * vsh * cosine) - + (susceptance * vsh * sine)); + const auto dpdAngle = + (conductance * voltage * vsh * sine) - (susceptance * voltage * vsh * cosine); + matrixDataValue.assignCheckCol(offset + ashIndex, busAngle, dpdAngle); + matrixDataValue.assign(offset + ashIndex, offset + ashIndex, -dpdAngle); + matrixDataValue.assign(offset + ashIndex, + offset + vshIndex, + (-conductance * voltage * cosine) - (susceptance * voltage * sine)); + matrixDataValue.assign(offset + ashIndex, offset + pshIndex, -1.0); + + matrixDataValue.assignCheckCol(offset + vshIndex, + busVoltage, + (-2.0 * susceptance * voltage) - (conductance * vsh * sine) + + (susceptance * vsh * cosine)); + const auto dqdAngle = + (-conductance * voltage * vsh * cosine) - (susceptance * voltage * vsh * sine); + matrixDataValue.assignCheckCol(offset + vshIndex, busAngle, dqdAngle); + matrixDataValue.assign(offset + vshIndex, offset + ashIndex, -dqdAngle); + matrixDataValue.assign(offset + vshIndex, + offset + vshIndex, + (-conductance * voltage * sine) + (susceptance * voltage * cosine)); + matrixDataValue.assign(offset + vshIndex, offset + qshIndex, -1.0); + + if ((control == Control::PQ) || (control == Control::PV)) { + matrixDataValue.assign(offset + pshIndex, offset + pshIndex, -1.0); + } else { + matrixDataValue.assignCheckCol(offset + pshIndex, dcVoltage, 1.0); + matrixDataValue.assignCheckCol(offset + pshIndex, dcReferenceVoltage, -1.0); + } + if ((control == Control::PQ) || (control == Control::VQ)) { + matrixDataValue.assign(offset + qshIndex, offset + qshIndex, -1.0); + } else { + matrixDataValue.assignCheckCol(offset + qshIndex, busVoltage, -1.0); + } + + matrixDataValue.assignCheckCol(offset + pdcIndex, + busVoltage, + (-conductance * vsh * cosine) + (susceptance * vsh * sine)); + const auto dpdcAngle = + (conductance * voltage * vsh * sine) + (susceptance * voltage * vsh * cosine); + matrixDataValue.assignCheckCol(offset + pdcIndex, busAngle, dpdcAngle); + matrixDataValue.assign(offset + pdcIndex, offset + ashIndex, -dpdcAngle); + matrixDataValue.assign(offset + pdcIndex, + offset + vshIndex, + (2.0 * conductance * vsh) - (conductance * voltage * cosine) + + (susceptance * voltage * sine)); + matrixDataValue.assign(offset + pdcIndex, offset + pdcIndex, 1.0); +} + +void VSCShunt::residual(const IOdata& inputs, + const StateData& stateDataValue, + double resid[], + const SolverMode& sMode) +{ + updateLocalCache(inputs, stateDataValue, sMode); + const auto offset = offsets.getAlgOffset(sMode); + const auto voltage = linkInfo.v1; + const auto angleDifference = linkInfo.theta1 - ash; + const auto cosine = std::cos(angleDifference); + const auto sine = std::sin(angleDifference); + const auto denominator = (r * r) + (x * x); + const auto conductance = r / denominator; + const auto susceptance = -x / denominator; + const auto dcDifference = + B2->getVoltage(stateDataValue, sMode) - dcReference->getVoltage(stateDataValue, sMode); + + resid[offset + ashIndex] = (conductance * voltage * voltage) - + (conductance * voltage * vsh * cosine) - (susceptance * voltage * vsh * sine) - psh; + resid[offset + vshIndex] = (-susceptance * voltage * voltage) - + (conductance * voltage * vsh * sine) + (susceptance * voltage * vsh * cosine) - qsh; + resid[offset + pshIndex] = + ((control == Control::PQ) || (control == Control::PV)) ? (p0 - psh) : (dcDifference - vdc0); + resid[offset + qshIndex] = + ((control == Control::PQ) || (control == Control::VQ)) ? (q0 - qsh) : (v0 - voltage); + resid[offset + pdcIndex] = (conductance * vsh * vsh) - (conductance * voltage * vsh * cosine) + + (susceptance * voltage * vsh * sine) + pdc; +} + +void VSCShunt::setState(CoreTime time, + const double state[], + const double /*dstateDt*/[], + const SolverMode& sMode) +{ + const auto offset = offsets.getAlgOffset(sMode); + ash = state[offset + ashIndex]; + vsh = state[offset + vshIndex]; + psh = state[offset + pshIndex]; + qsh = state[offset + qshIndex]; + pdc = state[offset + pdcIndex]; + prevTime = time; + updateLocalCache(); +} + +void VSCShunt::guessState(CoreTime /*time*/, + double state[], + double /*dstateDt*/[], + const SolverMode& sMode) +{ + const auto offset = offsets.getAlgOffset(sMode); + state[offset + ashIndex] = ash; + state[offset + vshIndex] = vsh; + state[offset + pshIndex] = psh; + state[offset + qshIndex] = qsh; + state[offset + pdcIndex] = pdc; +} + +void VSCShunt::getStateName(stringVec& stNames, + const SolverMode& sMode, + const std::string& prefix) const +{ + const auto offset = offsets.getAlgOffset(sMode); + const std::string statePrefix = prefix + getName() + ':'; + stNames[offset + ashIndex] = statePrefix + "ash"; + stNames[offset + vshIndex] = statePrefix + "vsh"; + stNames[offset + pshIndex] = statePrefix + "psh"; + stNames[offset + qshIndex] = statePrefix + "qsh"; + stNames[offset + pdcIndex] = statePrefix + "pdc"; +} + +double VSCShunt::getRealPower(id_type_t busId) const +{ + if ((busId == 3) || ((dcReference != nullptr) && (busId == dcReference->getID()))) { + const auto difference = dcVoltageDifference(); + if (std::abs(difference) < 1e-12) { + return 0.0; + } + if (currentBalance) { + return pdc / difference; + } + return pdc * static_cast(dcReference)->getVoltage() / difference; + } + return AcDcConverter::getRealPower(busId); +} + +double VSCShunt::getReactivePower(id_type_t busId) const +{ + if ((busId == 3) || ((dcReference != nullptr) && (busId == dcReference->getID()))) { + return 0.0; + } + return AcDcConverter::getReactivePower(busId); +} + +} // namespace griddyn::links diff --git a/src/griddyn/links/VSCShunt.h b/src/griddyn/links/VSCShunt.h new file mode 100644 index 000000000..4f2660a80 --- /dev/null +++ b/src/griddyn/links/VSCShunt.h @@ -0,0 +1,115 @@ +/* + * Copyright (c) 2014-2026, Lawrence Livermore National Security + * See the top-level NOTICE for additional details. All rights reserved. + * SPDX-License-Identifier: BSD-3-Clause + */ + +#pragma once + +#include "AcDcConverter.h" +#include +#include + +namespace griddyn { +class DcBus; + +namespace links { + /** + * @brief ANDES-compatible shunt voltage-source converter. + * + * This is a static AC/DC three-terminal link: terminal 1 is the AC bus and + * terminals 2 and 3 are the positive and negative DC nodes. Its algebraic + * equations follow ANDES' ``VSCShunt`` power-flow model. + */ + class VSCShunt final: public AcDcConverter { + public: + enum class Control { PQ = 0, PV = 1, VQ = 2, VV = 3 }; + + private: + DcBus* dcReference = nullptr; + Control control = Control::PQ; + double v0 = 1.0; + double p0 = 0.0; + double q0 = 0.0; + double vdc0 = 1.0; + double k0 = 0.0; + double k1 = 0.0; + double k2 = 0.0; + double droop = 0.0; + double droopK = 0.0; + double vhigh = 9999.0; + double vlow = 0.0; + double vshmax = 1.1; + double vshmin = 0.9; + double ishmax = 2.0; + bool currentBalance = false; + + // Saved algebraic variables: ash, vsh, psh, qsh, and pdc. + double ash = 0.0; + double vsh = 1.0; + double psh = 0.0; + double qsh = 0.0; + double pdc = 0.0; + + double dcVoltageDifference() const; + void updateFlows(double dcVoltage, double dcReferenceVoltage); + + public: + explicit VSCShunt(const std::string& objName = "vscshunt_$"); + + CoreObject* clone(CoreObject* obj = nullptr) const override; + void updateBus(GridBus* bus, index_t busnumber) override; + count_t terminalCount() const override { return 3; } + bool isConnected() const override; + GridBus* getBus(index_t busInd) const override; + void followNetwork(int network, std::queue& stk) override; + + void + set(std::string_view param, double val, units::unit unitType = units::defunit) override; + void set(std::string_view param, std::string_view val) override; + + void pFlowObjectInitializeA(CoreTime time0, std::uint32_t flags) override; + void dynObjectInitializeA(CoreTime time0, std::uint32_t flags) override; + StateSizes localStateSizes(const SolverMode& sMode) const override; + count_t localJacobianCount(const SolverMode& sMode) const override; + + void updateLocalCache() override; + void updateLocalCache(const IOdata& inputs, + const StateData& stateDataValue, + const SolverMode& sMode) override; + void ioPartialDerivatives(id_type_t busId, + const StateData& stateDataValue, + MatrixData& matrixDataValue, + const IOlocs& inputLocs, + const SolverMode& sMode) override; + void outputPartialDerivatives(id_type_t busId, + const StateData& stateDataValue, + MatrixData& matrixDataValue, + const SolverMode& sMode) override; + count_t outputDependencyCount(index_t num, const SolverMode& sMode) const override; + void jacobianElements(const IOdata& inputs, + const StateData& stateDataValue, + MatrixData& matrixDataValue, + const IOlocs& inputLocs, + const SolverMode& sMode) override; + void residual(const IOdata& inputs, + const StateData& stateDataValue, + double resid[], + const SolverMode& sMode) override; + void setState(CoreTime time, + const double state[], + const double dstateDt[], + const SolverMode& sMode) override; + void guessState(CoreTime time, + double state[], + double dstateDt[], + const SolverMode& sMode) override; + void getStateName(stringVec& stNames, + const SolverMode& sMode, + const std::string& prefix = "") const override; + double getRealPower(id_type_t busId = 0) const override; + double getReactivePower(id_type_t busId = 0) const override; + }; + +} // namespace links +} // namespace griddyn diff --git a/test/libraryTests/testJsonReader.cpp b/test/libraryTests/testJsonReader.cpp index f394c66e5..fd4bfbd11 100644 --- a/test/libraryTests/testJsonReader.cpp +++ b/test/libraryTests/testJsonReader.cpp @@ -7,20 +7,36 @@ // test case for CoreObject object #include "../gtestHelper.h" +#include "fileInput/fileInput.h" #include "formatInterpreters/jsonReaderElement.h" +#include "griddyn/links/DcLink.h" +#include "griddyn/links/VSCShunt.h" +#include "griddyn/primary/AcBus.h" +#include "griddyn/primary/DcBus.h" +#include +#include +#include #include #include +#include +#include #include #include static constexpr std::string_view elementReaderTestDirectory{GRIDDYN_TEST_DIRECTORY "/element_reader_tests/"}; +static constexpr std::string_view andesTestDirectory{GRIDDYN_TEST_DIRECTORY "/andes_tests/"}; static std::string makeElementReaderTestPath(std::string_view fileName) { return std::string{elementReaderTestDirectory} + std::string{fileName}; } +static std::string makeAndesTestPath(std::string_view fileName) +{ + return std::string{andesTestDirectory} + std::string{fileName}; +} + TEST(JsonElementReaderTests, JsonElementReaderTest1) { JsonReaderElement reader; @@ -193,3 +209,204 @@ TEST(JsonElementReaderTests, JsonElementReaderTest4) EXPECT_EQ(main->getName(), "main_element"); */ } + +TEST(AndesDcReaderTests, ImportsAllAndesDcComponents) +{ + auto simulation = std::make_unique(); + griddyn::loadFile(simulation.get(), makeAndesTestPath("andes_dc_components.json")); + + EXPECT_NE(dynamic_cast(simulation->find("ground_node")), nullptr); + EXPECT_NE(dynamic_cast(simulation->find("node1")), nullptr); + for (const auto* name : {"r", "l", "c", "rls", "rcp", "rlcp", "rcs", "rlcs"}) { + EXPECT_NE(dynamic_cast(simulation->find(name)), nullptr) << name; + } + auto* vsc = dynamic_cast(simulation->find("vsc")); + ASSERT_NE(vsc, nullptr); + EXPECT_EQ(vsc->terminalCount(), 3U); + ASSERT_NE(vsc->getBus(3), nullptr); + EXPECT_EQ(vsc->getBus(3)->getName(), "ground_node"); +} + +TEST(AndesVSCShuntTests, MatchesAndesPqReferencePoint) +{ + auto simulation = std::make_unique(); + auto* acBus = new griddyn::AcBus("ac"); + acBus->set("type", "swing"); + acBus->set("voltage", 1.0); + acBus->set("angle", 0.0); + auto* dcBus = new griddyn::DcBus("dc_positive"); + dcBus->set("type", "swing"); + dcBus->set("voltage", 1.0); + auto* dcReference = new griddyn::DcBus("dc_reference"); + dcReference->set("type", "swing"); + dcReference->set("voltage", 0.0); + auto* converter = new griddyn::links::VSCShunt("vsc"); + converter->set("r", 0.0025); + converter->set("x", 0.06); + converter->set("control", 0.0); // ANDES PQ mode + converter->set("p0", -0.1); + converter->set("q0", 0.0); + converter->set("andes_current_balance", 1.0); + converter->updateBus(acBus, 1); + converter->updateBus(dcBus, 2); + converter->updateBus(dcReference, 3); + simulation->add(acBus); + simulation->add(dcBus); + simulation->add(dcReference); + simulation->add(converter); + + EXPECT_EQ(simulation->powerflow(), 0); + // ANDES VSCShunt with rsh=0.0025, xsh=0.06, p0=-0.1 and q0=0 + // yields pdc=-0.100025 at Vdc=1. The DC terminals see -pdc/Vdc + // and pdc/Vdc, respectively. These reference values are evaluated + // from Andes' solved VSCShunt equations. + EXPECT_NEAR(converter->getRealPower(2), 0.100025, 1e-7); + EXPECT_NEAR(converter->getRealPower(3), -0.100025, 1e-7); + EXPECT_NEAR(converter->getRealPower(1), 0.1, 1e-9); + EXPECT_NEAR(converter->getReactivePower(1), 0.0, 1e-9); +} + +TEST(AndesVSCShuntTests, MatchesKundurVsc2OperatingPoint) +{ + // Boundary values and expected current were obtained from a one-worker + // Andes PFlow run of andes/cases/kundur/kundur_vsc.json. Pinning the + // external AC and DC voltages makes this a component-level comparison + // while the AC network importer is still being completed. + auto simulation = std::make_unique(); + auto* acBus = new griddyn::AcBus("ac"); + acBus->set("type", "swing"); + acBus->set("voltage", 0.968718231728735); + acBus->set("angle", 0.2893279182804385); + auto* dcBus = new griddyn::DcBus("dc_positive"); + dcBus->set("type", "swing"); + dcBus->set("voltage", 0.9989989369366359); + auto* dcReference = new griddyn::DcBus("dc_reference"); + dcReference->set("type", "swing"); + dcReference->set("voltage", 0.0); + auto* converter = new griddyn::links::VSCShunt("vsc"); + // Andes scales z=True parameters from VSC_2's 110 kV base to the + // 230 kV AC-bus base before evaluating its VSC equations. + converter->set("r", 0.0005718336483931947); + converter->set("x", 0.013724007561436672); + converter->set("control", 0.0); // Andes VSC_2, PQ mode + converter->set("p0", -0.1); + converter->set("q0", 0.0); + converter->set("andes_current_balance", 1.0); + converter->updateBus(acBus, 1); + converter->updateBus(dcBus, 2); + converter->updateBus(dcReference, 3); + simulation->add(acBus); + simulation->add(dcBus); + simulation->add(dcReference); + simulation->add(converter); + + ASSERT_EQ(simulation->powerflow(), 0); + // Andes: pdc=-0.10000609361072836, so -pdc / Vdc = 0.100106309... + EXPECT_NEAR(converter->getRealPower(2), 0.100106309, 1e-7); + EXPECT_NEAR(converter->getRealPower(3), -0.100106309, 1e-7); + EXPECT_NEAR(converter->getRealPower(1), 0.1, 1e-9); + EXPECT_NEAR(converter->getReactivePower(1), 0.0, 1e-9); +} + +TEST(AndesPowerFlowTests, MatchesCapturedKundurVscReference) +{ + std::ifstream input(makeAndesTestPath("andes_kundur_vsc_pflow_reference.json")); + ASSERT_TRUE(input.is_open()); + nlohmann::json reference; + input >> reference; + const auto tolerance = reference.at("tolerance").get(); + + auto simulation = std::make_unique(); + griddyn::loadFile(simulation.get(), makeAndesTestPath("andes_kundur_vsc_pflow.json")); + ASSERT_EQ(simulation->powerflow(), 0); + + for (std::size_t index = 0; index < reference.at("bus_voltage").size(); ++index) { + auto* bus = + dynamic_cast(simulation->find("Bus_" + std::to_string(index + 1))); + ASSERT_NE(bus, nullptr); + EXPECT_NEAR(bus->get("voltage"), reference["bus_voltage"][index].get(), tolerance); + EXPECT_NEAR(bus->get("angle"), reference["bus_angle"][index].get(), tolerance); + } + + for (std::size_t index = 0; index < reference.at("node_voltage").size(); ++index) { + auto* node = + dynamic_cast(simulation->find("Node_" + std::to_string(index))); + ASSERT_NE(node, nullptr); + EXPECT_NEAR(node->get("voltage"), + reference["node_voltage"][index].get(), + tolerance); + } + + for (index_t index = 0; index < 2; ++index) { + auto* converter = dynamic_cast( + simulation->find("VSC_" + std::to_string(index + 1))); + ASSERT_NE(converter, nullptr); + EXPECT_NEAR(converter->getRealPower(1), + -reference["vsc_psh"][index].get(), + tolerance); + EXPECT_NEAR(converter->getReactivePower(1), + -reference["vsc_qsh"][index].get(), + tolerance); + const auto dcVoltage = reference["node_voltage"][index + 1].get(); + EXPECT_NEAR(converter->getRealPower(2), + -reference["vsc_pdc"][index].get() / dcVoltage, + tolerance); + } +} + +TEST(AndesPowerFlowTests, MatchesCapturedTwoBusReference) +{ + std::ifstream input(makeAndesTestPath("andes_two_bus_pflow_reference.json")); + ASSERT_TRUE(input.is_open()); + nlohmann::json reference; + input >> reference; + const auto tolerance = reference.at("tolerance").get(); + + auto simulation = std::make_unique(); + griddyn::loadFile(simulation.get(), makeAndesTestPath("andes_two_bus_pflow.json")); + ASSERT_EQ(simulation->powerflow(), 0); + + const std::array busNames{"slack_bus", "load_bus"}; + for (std::size_t index = 0; index < busNames.size(); ++index) { + auto* bus = dynamic_cast(simulation->find(busNames[index])); + ASSERT_NE(bus, nullptr); + EXPECT_NEAR(bus->get("voltage"), reference["bus_voltage"][index].get(), tolerance); + EXPECT_NEAR(bus->get("angle"), reference["bus_angle"][index].get(), tolerance); + } +} + +TEST(AndesPowerFlowTests, MatchesCapturedVscResistorReference) +{ + std::ifstream input(makeAndesTestPath("andes_vsc_resistor_pflow_reference.json")); + ASSERT_TRUE(input.is_open()); + nlohmann::json reference; + input >> reference; + const auto tolerance = reference.at("tolerance").get(); + + auto simulation = std::make_unique(); + griddyn::loadFile(simulation.get(), makeAndesTestPath("andes_vsc_resistor_pflow.json")); + ASSERT_EQ(simulation->powerflow(), 0); + + auto* acBus = dynamic_cast(simulation->find("ac_slack")); + ASSERT_NE(acBus, nullptr); + EXPECT_NEAR(acBus->get("voltage"), reference["bus_voltage"][0].get(), tolerance); + EXPECT_NEAR(acBus->get("angle"), reference["bus_angle"][0].get(), tolerance); + + const std::array nodeNames{"dc_ground", "dc_node"}; + for (std::size_t index = 0; index < nodeNames.size(); ++index) { + auto* node = dynamic_cast(simulation->find(nodeNames[index])); + ASSERT_NE(node, nullptr); + EXPECT_NEAR(node->get("voltage"), + reference["node_voltage"][index].get(), + tolerance); + } + + auto* converter = dynamic_cast(simulation->find("vsc")); + ASSERT_NE(converter, nullptr); + EXPECT_NEAR(converter->getRealPower(1), -reference["vsc_psh"][0].get(), tolerance); + EXPECT_NEAR(converter->getReactivePower(1), -reference["vsc_qsh"][0].get(), tolerance); + const auto dcVoltage = reference["node_voltage"][1].get(); + EXPECT_NEAR(converter->getRealPower(2), + -reference["vsc_pdc"][0].get() / dcVoltage, + tolerance); +} diff --git a/test/test_files/andes_tests/andes_dc_components.json b/test/test_files/andes_tests/andes_dc_components.json new file mode 100644 index 000000000..efb673652 --- /dev/null +++ b/test/test_files/andes_tests/andes_dc_components.json @@ -0,0 +1,56 @@ +{ + "Bus": [{ "idx": 1, "name": "ac1", "Vn": 230.0, "v0": 1.0, "a0": 0.0 }], + "Node": [ + { "idx": 0, "name": "ground_node", "Vdcn": 100.0, "Idcn": 1.0, "v0": 0.0 }, + { "idx": 1, "name": "node1", "Vdcn": 100.0, "Idcn": 1.0, "v0": 1.0 }, + { "idx": 2, "name": "node2", "Vdcn": 100.0, "Idcn": 1.0, "v0": 1.0 } + ], + "Ground": [{ "idx": 0, "name": "ground", "node": 0, "voltage": 0.0 }], + "R": [{ "idx": 1, "name": "r", "node1": 1, "node2": 2, "R": 0.01 }], + "L": [{ "idx": 1, "name": "l", "node1": 1, "node2": 2, "L": 0.001 }], + "C": [{ "idx": 1, "name": "c", "node1": 1, "node2": 2, "C": 0.001 }], + "RLs": [ + { "idx": 1, "name": "rls", "node1": 1, "node2": 2, "R": 0.01, "L": 0.001 } + ], + "RCp": [ + { "idx": 1, "name": "rcp", "node1": 1, "node2": 2, "R": 0.01, "C": 0.001 } + ], + "RLCp": [ + { + "idx": 1, + "name": "rlcp", + "node1": 1, + "node2": 2, + "R": 0.01, + "L": 0.001, + "C": 0.001 + } + ], + "RCs": [ + { "idx": 1, "name": "rcs", "node1": 1, "node2": 2, "R": 0.01, "C": 0.001 } + ], + "RLCs": [ + { + "idx": 1, + "name": "rlcs", + "node1": 1, + "node2": 2, + "R": 0.01, + "L": 0.001, + "C": 0.001 + } + ], + "VSCShunt": [ + { + "idx": 1, + "name": "vsc", + "bus": 1, + "node1": 1, + "node2": 0, + "rsh": 0.0025, + "xsh": 0.06, + "p0": 0.1, + "Ishmax": 2.0 + } + ] +} diff --git a/test/test_files/andes_tests/andes_kundur_vsc_pflow.json b/test/test_files/andes_tests/andes_kundur_vsc_pflow.json new file mode 100644 index 000000000..7df1cab03 --- /dev/null +++ b/test/test_files/andes_tests/andes_kundur_vsc_pflow.json @@ -0,0 +1,183 @@ +{ + "Bus": [ + { "idx": 1, "Vn": 20.0, "v0": 1.0, "a0": 0.5702549171626017 }, + { "idx": 2, "Vn": 20.0, "v0": 0.99761, "a0": 0.3687461830443478 }, + { "idx": 3, "Vn": 20.0, "v0": 0.96263, "a0": 0.1853173146475028 }, + { "idx": 4, "Vn": 20.0, "v0": 0.81691, "a0": 0.462358662804314 }, + { "idx": 5, "Vn": 230.0, "v0": 0.97928, "a0": 0.4802029090767038 }, + { "idx": 6, "Vn": 230.0, "v0": 0.95796, "a0": 0.2838865294831329 }, + { "idx": 7, "Vn": 230.0, "v0": 0.9362, "a0": 0.1269011445832536 }, + { "idx": 8, "Vn": 230.0, "v0": 0.87904, "a0": -0.0805923235400888 }, + { "idx": 9, "Vn": 230.0, "v0": 0.89054, "a0": 0.09361771574772226 }, + { "idx": 10, "Vn": 230.0, "v0": 0.82958, "a0": 0.3366007088811167 } + ], + "PQ": [ + { "idx": "PQ_0", "bus": 7, "p0": 11.59, "q0": -0.735 }, + { "idx": "PQ_1", "bus": 8, "p0": 15.75, "q0": -0.899 } + ], + "PV": [ + { "idx": 2, "bus": 2, "p0": 7.0, "q0": 3.0, "v0": 1.0 }, + { "idx": 3, "bus": 3, "p0": 7.0, "q0": 5.5, "v0": 1.0 }, + { "idx": 4, "bus": 4, "p0": 7.0, "q0": -1.0, "v0": 1.0 } + ], + "Slack": [ + { + "idx": 1, + "bus": 1, + "p0": 7.45861, + "q0": 1.43612, + "v0": 1.0, + "a0": 0.5702549171626017 + } + ], + "Line": [ + { + "idx": "Line_0", + "bus1": 5, + "bus2": 6, + "r": 0.005, + "x": 0.05, + "b": 0.075 + }, + { + "idx": "Line_1", + "bus1": 5, + "bus2": 6, + "r": 0.00501, + "x": 0.05001, + "b": 0.075 + }, + { "idx": "Line_2", "bus1": 6, "bus2": 7, "r": 0.002, "x": 0.02, "b": 0.03 }, + { + "idx": "Line_3", + "bus1": 6, + "bus2": 7, + "r": 0.00201, + "x": 0.02001, + "b": 0.03 + }, + { + "idx": "Line_4", + "bus1": 7, + "bus2": 8, + "r": 0.02201, + "x": 0.22001, + "b": 0.33 + }, + { + "idx": "Line_5", + "bus1": 7, + "bus2": 8, + "r": 0.02202, + "x": 0.22002, + "b": 0.33 + }, + { "idx": "Line_6", "bus1": 7, "bus2": 8, "r": 0.022, "x": 0.22, "b": 0.33 }, + { "idx": "Line_7", "bus1": 8, "bus2": 9, "r": 0.002, "x": 0.02, "b": 0.03 }, + { + "idx": "Line_8", + "bus1": 8, + "bus2": 9, + "r": 0.00201, + "x": 0.02001, + "b": 0.03 + }, + { + "idx": "Line_9", + "bus1": 9, + "bus2": 10, + "r": 0.005, + "x": 0.05, + "b": 0.075 + }, + { + "idx": "Line_10", + "bus1": 9, + "bus2": 10, + "r": 0.00501, + "x": 0.05001, + "b": 0.075 + }, + { + "idx": "Line_11", + "bus1": 1, + "bus2": 5, + "r": 0.001, + "x": 0.012, + "b": 0.0 + }, + { + "idx": "Line_12", + "bus1": 2, + "bus2": 6, + "r": 0.001, + "x": 0.012, + "b": 0.0 + }, + { + "idx": "Line_13", + "bus1": 3, + "bus2": 9, + "r": 0.001, + "x": 0.012, + "b": 0.0 + }, + { + "idx": "Line_14", + "bus1": 4, + "bus2": 10, + "r": 0.001, + "x": 0.012, + "b": 0.0 + } + ], + "Node": [ + { "idx": 0, "name": "Node_0", "Vdcn": 100.0, "v0": 0.0 }, + { "idx": 1, "name": "Node_1", "Vdcn": 100.0, "v0": 1.0 }, + { "idx": 2, "name": "Node_2", "Vdcn": 100.0, "v0": 1.0 } + ], + "Ground": [{ "idx": 0, "node": 0, "voltage": 0.0 }], + "R": [{ "idx": 1, "name": "R1_2", "node1": 1, "node2": 2, "R": 0.01 }], + "VSCShunt": [ + { + "idx": 1, + "name": "VSC_1", + "bus": 1, + "node1": 1, + "node2": 0, + "Vn": 110.0, + "rsh": 0.0025, + "xsh": 0.06, + "control": 2.0, + "v0": 1.0, + "p0": 0.0, + "q0": 0.01, + "vdc0": 1.0, + "k0": 0.0, + "k1": 0.0, + "k2": 0.0, + "droop": 0.0, + "K": 0.0 + }, + { + "idx": 2, + "name": "VSC_2", + "bus": 6, + "node1": 2, + "node2": 0, + "Vn": 110.0, + "rsh": 0.0025, + "xsh": 0.06, + "control": 0.0, + "v0": 1.0, + "p0": -0.1, + "q0": 0.0, + "vdc0": 1.0, + "k0": 0.0, + "k1": 0.0, + "k2": 0.0, + "droop": 0.0, + "K": 0.0 + } + ] +} diff --git a/test/test_files/andes_tests/andes_kundur_vsc_pflow_reference.json b/test/test_files/andes_tests/andes_kundur_vsc_pflow_reference.json new file mode 100644 index 000000000..9d98c2912 --- /dev/null +++ b/test/test_files/andes_tests/andes_kundur_vsc_pflow_reference.json @@ -0,0 +1,20 @@ +{ + "source": "ANDES pflow: cases/kundur/kundur_vsc.json (single worker)", + "tolerance": 1e-6, + "bus_voltage": [ + 1.0000000000182983, 1.0000000000359046, 1.000000000019338, + 1.000000000008106, 0.9830217891686318, 0.968718231728735, + 0.9558683192567611, 0.9539200692156563, 0.9685248386762372, + 0.9837586556613533 + ], + "bus_angle": [ + 0.5702549171702109, 0.37375511439598824, 0.19134501468615742, + 0.37330102249881925, 0.48126395881854245, 0.2893279182804385, + 0.13822340710250802, -0.04156595505786644, 0.10691756193816289, + 0.28889347555716843 + ], + "node_voltage": [0.0, 1.0, 0.9989989369366359], + "vsc_psh": [0.1008835414369301, -0.1], + "vsc_qsh": [0.01, 0.0], + "vsc_pdc": [0.10010630633641154, -0.10000609361072836] +} diff --git a/test/test_files/andes_tests/andes_two_bus_pflow.json b/test/test_files/andes_tests/andes_two_bus_pflow.json new file mode 100644 index 000000000..797cf2bc3 --- /dev/null +++ b/test/test_files/andes_tests/andes_two_bus_pflow.json @@ -0,0 +1,21 @@ +{ + "Bus": [ + { "idx": 1, "name": "slack_bus", "Vn": 230.0, "v0": 1.0, "a0": 0.0 }, + { "idx": 2, "name": "load_bus", "Vn": 230.0, "v0": 0.98, "a0": -0.1 } + ], + "Slack": [{ "idx": 1, "bus": 1, "p0": 0.0, "q0": 0.0, "v0": 1.0, "a0": 0.0 }], + "PQ": [{ "idx": "PQ_1", "bus": 2, "p0": 1.0, "q0": 0.3 }], + "Line": [ + { + "idx": "Line_1", + "bus1": 1, + "bus2": 2, + "Vn1": 230.0, + "Vn2": 230.0, + "r": 0.01, + "x": 0.1, + "b": 0.02 + } + ], + "Node": [] +} diff --git a/test/test_files/andes_tests/andes_two_bus_pflow_reference.json b/test/test_files/andes_tests/andes_two_bus_pflow_reference.json new file mode 100644 index 000000000..3ed93dc1c --- /dev/null +++ b/test/test_files/andes_tests/andes_two_bus_pflow_reference.json @@ -0,0 +1,6 @@ +{ + "source": "ANDES pflow: andes_two_bus_pflow.json (single worker)", + "tolerance": 1e-6, + "bus_voltage": [1.0000000000000588, 0.9538230766211548], + "bus_angle": [5.872548199249326e-15, -0.10196800834019634] +} diff --git a/test/test_files/andes_tests/andes_vsc_resistor_pflow.json b/test/test_files/andes_tests/andes_vsc_resistor_pflow.json new file mode 100644 index 000000000..c9a6ad4f8 --- /dev/null +++ b/test/test_files/andes_tests/andes_vsc_resistor_pflow.json @@ -0,0 +1,27 @@ +{ + "Bus": [{ "idx": 1, "name": "ac_slack", "Vn": 230.0, "v0": 1.0, "a0": 0.0 }], + "Slack": [{ "idx": 1, "bus": 1, "p0": 0.0, "q0": 0.0, "v0": 1.0, "a0": 0.0 }], + "Node": [ + { "idx": 0, "name": "dc_ground", "Vdcn": 100.0, "v0": 0.0 }, + { "idx": 1, "name": "dc_node", "Vdcn": 100.0, "v0": 1.0 } + ], + "Ground": [{ "idx": 0, "node": 0, "voltage": 0.0 }], + "R": [{ "idx": 1, "name": "dc_resistor", "node1": 1, "node2": 0, "R": 10.0 }], + "VSCShunt": [ + { + "idx": 1, + "name": "vsc", + "bus": 1, + "node1": 1, + "node2": 0, + "Vn": 230.0, + "rsh": 0.0025, + "xsh": 0.06, + "control": 0.0, + "v0": 1.0, + "p0": 0.1, + "q0": 0.02, + "vdc0": 1.0 + } + ] +} diff --git a/test/test_files/andes_tests/andes_vsc_resistor_pflow_reference.json b/test/test_files/andes_tests/andes_vsc_resistor_pflow_reference.json new file mode 100644 index 000000000..f6c5578d6 --- /dev/null +++ b/test/test_files/andes_tests/andes_vsc_resistor_pflow_reference.json @@ -0,0 +1,10 @@ +{ + "source": "ANDES pflow: andes_vsc_resistor_pflow.json (single worker)", + "tolerance": 1e-6, + "bus_voltage": [1.0], + "bus_angle": [0.0], + "node_voltage": [1.2999924473942318e-20, 0.9998700013499972], + "vsc_psh": [0.1], + "vsc_qsh": [0.02], + "vsc_pdc": [0.0999740002686917] +}