diff --git a/CHANGELOG.md b/CHANGELOG.md
index 744f6ec98..7e3336671 100644
--- a/CHANGELOG.md
+++ b/CHANGELOG.md
@@ -82,6 +82,7 @@
 - Added `REECB` electrical-control model implementation for PhasorDynamics.
 - Added cases and validation for ACTIVSg10k, ACTIVSg200, ACTIVSg500, and WECC240.
 - Added IDA option to choose the consistent initial condition calculation type.
+- Implemented `tagDifferentiable()` for `PowerElectronics` models.
 
 ## v0.1
 
diff --git a/GridKit/Model/PowerElectronics/Capacitor/Capacitor.cpp b/GridKit/Model/PowerElectronics/Capacitor/Capacitor.cpp
index aab348679..6c0b5fdf6 100644
--- a/GridKit/Model/PowerElectronics/Capacitor/Capacitor.cpp
+++ b/GridKit/Model/PowerElectronics/Capacitor/Capacitor.cpp
@@ -49,6 +49,8 @@ namespace GridKit
   template <class ScalarT, typename IdxT>
   int Capacitor<ScalarT, IdxT>::tagDifferentiable()
   {
+    // All variables are differentials
+    std::fill(tag_.begin(), tag_.end(), true);
     return 0;
   }
 
diff --git a/GridKit/Model/PowerElectronics/CircuitComponent.hpp b/GridKit/Model/PowerElectronics/CircuitComponent.hpp
index c6fcad6dd..62e67dc84 100644
--- a/GridKit/Model/PowerElectronics/CircuitComponent.hpp
+++ b/GridKit/Model/PowerElectronics/CircuitComponent.hpp
@@ -132,10 +132,12 @@ namespace GridKit
       jacobian_coo_cols_   = std::make_unique<IdxT[]>(static_cast<size_t>(nnz_));
       jacobian_coo_values_ = std::make_unique<RealT[]>(static_cast<size_t>(nnz_));
 
-      y_ext_            = std::make_unique<const ScalarT*[]>(static_cast<size_t>(size_));
-      yp_ext_           = std::make_unique<const ScalarT*[]>(static_cast<size_t>(size_));
-      f_ext_            = std::make_unique<ScalarT*[]>(static_cast<size_t>(size_));
-      connection_nodes_ = std::make_unique<IdxT[]>(static_cast<size_t>(size_));
+      y_ext_            = std::make_unique<const ScalarT*[]>(static_cast<IdxT>(size_));
+      yp_ext_           = std::make_unique<const ScalarT*[]>(static_cast<IdxT>(size_));
+      f_ext_            = std::make_unique<ScalarT*[]>(static_cast<IdxT>(size_));
+      connection_nodes_ = std::make_unique<IdxT[]>(static_cast<IdxT>(size_));
+
+      tag_.resize(static_cast<IdxT>(size_));
 
       if (!allocated_)
       {
diff --git a/GridKit/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.cpp b/GridKit/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.cpp
index dfa6f1a98..aa538933f 100644
--- a/GridKit/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.cpp
+++ b/GridKit/Model/PowerElectronics/DistributedGenerator/DistributedGenerator.cpp
@@ -3,7 +3,6 @@
 #include "DistributedGenerator.hpp"
 
 #include <cmath>
-#include <iostream>
 #include <vector>
 
 namespace GridKit
@@ -74,6 +73,8 @@ namespace GridKit
   template <class ScalarT, typename IdxT>
   int DistributedGenerator<ScalarT, IdxT>::tagDifferentiable()
   {
+    // All variables are differentials
+    std::fill(tag_.begin(), tag_.end(), true);
     return 0;
   }
 
diff --git a/GridKit/Model/PowerElectronics/InductionMotor/InductionMotor.cpp b/GridKit/Model/PowerElectronics/InductionMotor/InductionMotor.cpp
index 0224e6ebe..74aabb549 100644
--- a/GridKit/Model/PowerElectronics/InductionMotor/InductionMotor.cpp
+++ b/GridKit/Model/PowerElectronics/InductionMotor/InductionMotor.cpp
@@ -60,6 +60,8 @@ namespace GridKit
   template <class ScalarT, typename IdxT>
   int InductionMotor<ScalarT, IdxT>::tagDifferentiable()
   {
+    // All variables are differentials
+    std::fill(tag_.begin(), tag_.end(), true);
     return 0;
   }
 
diff --git a/GridKit/Model/PowerElectronics/Inductor/Inductor.cpp b/GridKit/Model/PowerElectronics/Inductor/Inductor.cpp
index b192cebc5..0e6e7ea28 100644
--- a/GridKit/Model/PowerElectronics/Inductor/Inductor.cpp
+++ b/GridKit/Model/PowerElectronics/Inductor/Inductor.cpp
@@ -48,6 +48,8 @@ namespace GridKit
   template <class ScalarT, typename IdxT>
   int Inductor<ScalarT, IdxT>::tagDifferentiable()
   {
+    // All variables are differentials
+    std::fill(tag_.begin(), tag_.end(), true);
     return 0;
   }
 
diff --git a/GridKit/Model/PowerElectronics/LinearTransformer/LinearTransformer.cpp b/GridKit/Model/PowerElectronics/LinearTransformer/LinearTransformer.cpp
index fcfc4a02d..455ef202f 100644
--- a/GridKit/Model/PowerElectronics/LinearTransformer/LinearTransformer.cpp
+++ b/GridKit/Model/PowerElectronics/LinearTransformer/LinearTransformer.cpp
@@ -62,6 +62,8 @@ namespace GridKit
   template <class ScalarT, typename IdxT>
   int LinearTransformer<ScalarT, IdxT>::tagDifferentiable()
   {
+    // All variables are differentials
+    std::fill(tag_.begin(), tag_.end(), true);
     return 0;
   }
 
diff --git a/GridKit/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.cpp b/GridKit/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.cpp
index ffc02dfb3..462c200f8 100644
--- a/GridKit/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.cpp
+++ b/GridKit/Model/PowerElectronics/MicrogridBusDQ/MicrogridBusDQ.cpp
@@ -47,9 +47,7 @@ namespace GridKit
     return 0;
   }
 
-  /*
-   * \brief Identify differential variables
-   */
+  /// There are no internal variables in this component, so \ref tag_ can be set arbitrarily.
   template <class ScalarT, typename IdxT>
   int MicrogridBusDQ<ScalarT, IdxT>::tagDifferentiable()
   {
diff --git a/GridKit/Model/PowerElectronics/MicrogridLine/MicrogridLine.cpp b/GridKit/Model/PowerElectronics/MicrogridLine/MicrogridLine.cpp
index 65c3fb0a0..efa5c2a36 100644
--- a/GridKit/Model/PowerElectronics/MicrogridLine/MicrogridLine.cpp
+++ b/GridKit/Model/PowerElectronics/MicrogridLine/MicrogridLine.cpp
@@ -62,6 +62,8 @@ namespace GridKit
   template <class ScalarT, typename IdxT>
   int MicrogridLine<ScalarT, IdxT>::tagDifferentiable()
   {
+    // All variables are differentials
+    std::fill(tag_.begin(), tag_.end(), true);
     return 0;
   }
 
diff --git a/GridKit/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.cpp b/GridKit/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.cpp
index a14a038f1..36ac6c375 100644
--- a/GridKit/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.cpp
+++ b/GridKit/Model/PowerElectronics/MicrogridLoad/MicrogridLoad.cpp
@@ -58,6 +58,8 @@ namespace GridKit
   template <class ScalarT, typename IdxT>
   int MicrogridLoad<ScalarT, IdxT>::tagDifferentiable()
   {
+    // All variables are differentials
+    std::fill(tag_.begin(), tag_.end(), true);
     return 0;
   }
 
diff --git a/GridKit/Model/PowerElectronics/NodeBase.hpp b/GridKit/Model/PowerElectronics/NodeBase.hpp
index 3c7f7c80b..077854b97 100644
--- a/GridKit/Model/PowerElectronics/NodeBase.hpp
+++ b/GridKit/Model/PowerElectronics/NodeBase.hpp
@@ -145,7 +145,7 @@ namespace GridKit
           allocateVectors(static_cast<IdxT>(size));
         }
 
-        tag_.resize(size);
+        tag_.resize(size, false);
         variable_indices_.resize(size);
         residual_indices_.resize(size);
 
diff --git a/GridKit/Model/PowerElectronics/Resistor/Resistor.cpp b/GridKit/Model/PowerElectronics/Resistor/Resistor.cpp
index 5c6c4a1cb..8fba09a55 100644
--- a/GridKit/Model/PowerElectronics/Resistor/Resistor.cpp
+++ b/GridKit/Model/PowerElectronics/Resistor/Resistor.cpp
@@ -43,9 +43,7 @@ namespace GridKit
     return 0;
   }
 
-  /*
-   * \brief Identify differential variables
-   */
+  /// There are no internal variables in this component, so \ref tag_ can be set arbitrarily.
   template <class ScalarT, typename IdxT>
   int Resistor<ScalarT, IdxT>::tagDifferentiable()
   {
diff --git a/GridKit/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.cpp b/GridKit/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.cpp
index c3a647762..5ec5b3b91 100644
--- a/GridKit/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.cpp
+++ b/GridKit/Model/PowerElectronics/SynchronousMachine/SynchronousMachine.cpp
@@ -77,6 +77,8 @@ namespace GridKit
   template <class ScalarT, typename IdxT>
   int SynchronousMachine<ScalarT, IdxT>::tagDifferentiable()
   {
+    // All variables are differentials
+    std::fill(tag_.begin(), tag_.end(), true);
     return 0;
   }
 
diff --git a/GridKit/Model/PowerElectronics/SystemModelPowerElectronics.hpp b/GridKit/Model/PowerElectronics/SystemModelPowerElectronics.hpp
index 57b705bbb..51c8dcee1 100644
--- a/GridKit/Model/PowerElectronics/SystemModelPowerElectronics.hpp
+++ b/GridKit/Model/PowerElectronics/SystemModelPowerElectronics.hpp
@@ -292,8 +292,55 @@ namespace GridKit
       return CircuitComponent<ScalarT, IdxT>::initialize();
     }
 
+    /**
+     * @brief Tags all system variables as differentiable, based on what the
+     * components that own those variables tag them as.
+     *
+     * Starts by asking all components to tag their differentiables. This implementation
+     * assumes all node variables are algebraic, and will not ask nodes to tag their differentiables.
+     * Sets all variables to algebraic (`false`) to start, then loops over all component internal variables.
+     * Re-creates the same internal variable to system variables mapping as in \ref allocate() - all
+     * internal variables from the same component are stored contiguously in a block, and blocks are
+     * stored contiguously in the same order as \ref components_, with node variables at the end.
+     * Each internal variable's tag in the system is set to its tag in the component.
+     */
     int tagDifferentiable() final
     {
+      // Ask all component to tag their differentiables
+      for (size_t i = 0; i < components_.size(); i++)
+      {
+        component_type* component = components_[i];
+
+        // Bubble up errors if necessary
+        if (int err = component->tagDifferentiable())
+        {
+          return err;
+        }
+      }
+
+      // Fill tags with a default value (false) for node variables. Assumed to be algebraic here.
+      std::fill(tag_.begin(), tag_.end(), false);
+
+      // Copy tags for internal variables from their components - going in the order as described above
+      size_t idx = 0;
+      for (component_type* comp : components_)
+      {
+        const auto& external_indices = comp->getExternIndices();
+
+        // Loop over all component variables - including externals
+        for (IdxT i = 0; i < comp->size(); i++)
+        {
+          // Discard externals
+          if (!external_indices.contains(i))
+          {
+            tag_[idx] = comp->tag()[i];
+
+            // Ensures internal variables are contiguous, and in the same order as the component
+            idx++;
+          }
+        }
+      }
+
       return 0;
     }
 
diff --git a/GridKit/Model/PowerElectronics/TransmissionLine/TransmissionLine.cpp b/GridKit/Model/PowerElectronics/TransmissionLine/TransmissionLine.cpp
index 3adbbef58..f4700e275 100644
--- a/GridKit/Model/PowerElectronics/TransmissionLine/TransmissionLine.cpp
+++ b/GridKit/Model/PowerElectronics/TransmissionLine/TransmissionLine.cpp
@@ -60,6 +60,8 @@ namespace GridKit
   template <class ScalarT, typename IdxT>
   int TransmissionLine<ScalarT, IdxT>::tagDifferentiable()
   {
+    // All variables are algebraics
+    std::fill(tag_.begin(), tag_.end(), false);
     return 0;
   }
 
diff --git a/GridKit/Model/PowerElectronics/VoltageSource/VoltageSource.cpp b/GridKit/Model/PowerElectronics/VoltageSource/VoltageSource.cpp
index 9a4accd0b..52f021a5a 100644
--- a/GridKit/Model/PowerElectronics/VoltageSource/VoltageSource.cpp
+++ b/GridKit/Model/PowerElectronics/VoltageSource/VoltageSource.cpp
@@ -49,6 +49,8 @@ namespace GridKit
   template <class ScalarT, typename IdxT>
   int VoltageSource<ScalarT, IdxT>::tagDifferentiable()
   {
+    // All variables are algebraics
+    std::fill(tag_.begin(), tag_.end(), false);
     return 0;
   }
 
diff --git a/examples/PowerElectronics/Microgrid/Microgrid.cpp b/examples/PowerElectronics/Microgrid/Microgrid.cpp
index 74497fc41..dfd71a7e1 100644
--- a/examples/PowerElectronics/Microgrid/Microgrid.cpp
+++ b/examples/PowerElectronics/Microgrid/Microgrid.cpp
@@ -223,6 +223,33 @@ int main(int /* argc */, char const** /* argv */)
     sysmodel->getCsrJacobian()->print();
   }
 
+  sysmodel->tagDifferentiable();
+
+  bool all_internal_diff = true;
+  bool all_external_alg  = true;
+
+  const size_t num_node_vars = bus1.size() + bus2.size() + bus3.size() + bus4.size() + dg_signal.size();
+
+  for (size_t i = 0; i < sysmodel->size() - num_node_vars; i++)
+  {
+    all_internal_diff = all_internal_diff && sysmodel->tag()[i];
+    if (!sysmodel->tag()[i])
+    {
+      std::cout << "Unexepected algebraic-tagged internal variable found in index " << i << '\n';
+    }
+  }
+
+  for (size_t i = sysmodel->size() - num_node_vars; i < sysmodel->size(); i++)
+  {
+    all_external_alg = all_external_alg && !sysmodel->tag()[i];
+    if (sysmodel->tag()[i])
+    {
+      std::cout << "Unexepected differential-tagged external variable found in index " << i << '\n';
+    }
+  }
+
+  std::cout << "Verify all internal variables are differential: " << all_internal_diff << ", and all external variables are algebraic: " << all_external_alg << '\n';
+
   // Create numerical integrator and configure it for the generator model
   auto* idas = new AnalysisManager::Sundials::Ida<double, size_t>(sysmodel);