diff --git a/src/include/serif/solver/seed/fitted.hpp b/src/include/serif/solver/seed/fitted.hpp
new file mode 100644
index 00000000..5a8872c0
--- /dev/null
+++ b/src/include/serif/solver/seed/fitted.hpp
@@ -0,0 +1,50 @@
+#pragma once
+
+#include "serif/solver/seed/polytrope.hpp"
+
+namespace serif::solver::seed {
+    /**
+     * @brief An explicitly stretched initial guess, distinct from an equilibrium Lane–Emden scaling.
+     * @note The source EOS, density and enthalpy stay unchanged at corresponding fractional radii.
+     * Gravity follows the stretched density; this generally breaks hydrostatic balance, intentionally.
+     */
+    struct FittedRadialSeed final {
+        RadialSeed source;
+        dimensions::LengthValue radius;
+        dimensions::MassValue mass;
+        double moment_of_inertia{};
+        double angular_velocity{};
+    };
+
+    /** @brief Evaluate the original equilibrium profile without changing its physical scaling. */
+    [[nodiscard]] inline auto evaluate_radial_seed(const RadialSeed& seed, dimensions::LengthValue radius) {
+        return utils::laneemden::evaluate_physical(seed.solution, seed.scales, radius);
+    }
+
+    /**
+     * @brief Copy a radial seed onto another radius without changing its material thermodynamics.
+     * @note For stretch s, enclosed mass scales as s^3, potential as s^2, gravity as s and inertia as s^5.
+     * Omega becomes Omega_source/s^5, preserving the angular momentum inferred from the source guess.
+     * The original PhysicalScales remain intact; no equilibrium or target-constraint satisfaction is implied.
+     */
+    [[nodiscard]] std::expected<FittedRadialSeed, std::string> fit_radial_seed(const RadialSeed& seed, dimensions::LengthValue radius);
+
+    /** @brief Evaluate the fitted density and its consistent spherical gravity, retaining the source enthalpy. */
+    [[nodiscard]] std::expected<utils::laneemden::PhysicalSample, std::string> evaluate_radial_seed(
+        const FittedRadialSeed& seed, dimensions::LengthValue radius);
+
+    namespace detail {
+        /** @brief Inspect the unchanged physical source of a strict seed. */
+        [[nodiscard]] inline const RadialSeed& radial_source(const RadialSeed& seed) { return seed; }
+        /** @brief Inspect the unchanged physical source of an explicitly fitted guess. */
+        [[nodiscard]] inline const RadialSeed& radial_source(const FittedRadialSeed& seed) { return seed.source; }
+        /** @brief Use the equilibrium radius for ordinary projection. */
+        [[nodiscard]] inline dimensions::LengthValue radial_radius(const RadialSeed& seed) { return seed.scales.radius; }
+        /** @brief Use the explicit destination radius for fitted projection. */
+        [[nodiscard]] inline dimensions::LengthValue radial_radius(const FittedRadialSeed& seed) { return seed.radius; }
+        /** @brief Retain the original spherical mass for ordinary projection. */
+        [[nodiscard]] inline dimensions::MassValue radial_mass(const RadialSeed& seed) { return seed.scales.mass; }
+        /** @brief Use the mass of the stretched density for fitted gravity projection. */
+        [[nodiscard]] inline dimensions::MassValue radial_mass(const FittedRadialSeed& seed) { return seed.mass; }
+    }
+}
diff --git a/src/include/serif/solver/seed/polytrope.hpp b/src/include/serif/solver/seed/polytrope.hpp
new file mode 100644
index 00000000..122fb79c
--- /dev/null
+++ b/src/include/serif/solver/seed/polytrope.hpp
@@ -0,0 +1,54 @@
+#pragma once
+
+#include <array>
+#include <concepts>
+#include <expected>
+#include <optional>
+#include <string>
+
+#include "serif/eos/models/polytropic.hpp"
+#include "serif/model/stellar_model.hpp"
+#include "serif/model/specifications/fixed_total_mass.hpp"
+#include "serif/model/specifications/fixed_angular_momentum.hpp"
+#include "serif/utils/laneemden/laneemden.hpp"
+
+namespace serif::solver::seed {
+    /** @brief Numerical choices for an EOS-consistent spherical starting guess, separate from solver normalization. */
+    struct PolytropeOptions final {
+        dimensions::QuantityValue<dimensions::GravitationalConstant> gravitational_constant{1.0};
+        std::optional<dimensions::DensityValue> central_density;
+        utils::laneemden::IntegrationOptions integration;
+    };
+
+    /** @brief Own a spherical physical profile and its first rotation estimate without changing the stellar model. */
+    struct RadialSeed final {
+        utils::laneemden::Solution solution;
+        utils::laneemden::PhysicalScales scales;
+        double moment_of_inertia{};
+        double angular_velocity{};
+        std::array<double, 3> center{};
+    };
+
+    /**
+     * @brief Fix the EOS K and requested mass, solving for central density when n differs from three.
+     * @note At n=3 mass does not determine central density: supply that homology scale explicitly.
+     * An explicit density specifies the initial guess only; its mass need not satisfy the final constraint.
+     */
+    [[nodiscard]] std::expected<RadialSeed, std::string> make_polytropic_seed(
+        const eos::models::Polytrope& equation_of_state, dimensions::MassValue mass, PolytropeOptions options = {});
+
+    /** @brief Read the authoritative EOS, mass and optional rotation prescription from a model. */
+    template <model::StellarModelType Model>
+        requires std::same_as<typename Model::EquationOfState, eos::models::Polytrope> && Model::template contains_specification<model::FixedTotalMass>
+    [[nodiscard]] std::expected<RadialSeed, std::string> make_polytropic_seed(const Model& model, PolytropeOptions options = {}) {
+        auto seed = make_polytropic_seed(model.equation_of_state(), model.template specification<model::FixedTotalMass>().target(), options);
+        if (seed) {
+            if constexpr (Model::template contains_specification<model::FixedAngularMomentum>) {
+                const auto& rotation = model.template specification<model::FixedAngularMomentum>();
+                seed->angular_velocity = rotation.target().value() / seed->moment_of_inertia;
+                seed->center = rotation.center();
+            }
+        }
+        return seed;
+    }
+}
diff --git a/src/include/serif/solver/seed/state.hpp b/src/include/serif/solver/seed/state.hpp
new file mode 100644
index 00000000..79e19b5d
--- /dev/null
+++ b/src/include/serif/solver/seed/state.hpp
@@ -0,0 +1,9 @@
+#pragma once
+
+#include "serif/backend/mfem/solver/seed/state.hpp"
+#include "serif/backend/mfem/solver/seed/fitted.hpp"
+
+namespace serif::solver::seed {
+    using backend::mfem::make_polytropic_state;
+    using backend::mfem::fit_radial_seed;
+}
diff --git a/src/lib/serif/solver/seed/fitted.cpp b/src/lib/serif/solver/seed/fitted.cpp
new file mode 100644
index 00000000..b9fce5a1
--- /dev/null
+++ b/src/lib/serif/solver/seed/fitted.cpp
@@ -0,0 +1,87 @@
+#include "serif/solver/seed/fitted.hpp"
+
+#include <cmath>
+
+namespace serif::solver::seed {
+    /** @brief Validate a regular source and own its explicit, thermodynamically unchanged radial stretch. */
+    std::expected<FittedRadialSeed, std::string> fit_radial_seed(const RadialSeed& seed, dimensions::LengthValue radius) {
+        const auto& samples = seed.solution.samples;
+        if (samples.size() < 2 || samples.front().xi != 0.0 || samples.front().theta != 1.0 || samples.front().derivative != 0.0 ||
+            !seed.solution.surface || samples.back().xi != seed.solution.surface->xi || samples.back().theta != 0.0 ||
+            samples.back().derivative != seed.solution.surface->derivative) {
+            return std::unexpected("Fitting a radial guess requires a complete regular Lane–Emden profile ending at its finite surface.");
+        }
+        double previous_radius = -1.0;
+        for (const auto& sample : samples) {
+            if (!std::isfinite(sample.xi) || sample.xi <= previous_radius || !std::isfinite(sample.theta) || sample.theta < 0.0 || !std::isfinite(sample.derivative)) {
+                return std::unexpected("Fitting a radial guess requires finite nonnegative samples at strictly increasing dimensionless radii.");
+            }
+            previous_radius = sample.xi;
+        }
+        const auto center = evaluate_radial_seed(seed, dimensions::LengthValue{0.0});
+        if (!center) {
+            return std::unexpected(center.error());
+        }
+        if (!std::isfinite(radius.value()) || radius.value() <= 0.0 || !std::isfinite(seed.moment_of_inertia) ||
+            seed.moment_of_inertia <= 0.0 || !std::isfinite(seed.angular_velocity)) {
+            return std::unexpected("Fitting a radial guess requires a finite positive radius and inertia and a finite angular velocity.");
+        }
+        for (const double coordinate : seed.center) {
+            if (!std::isfinite(coordinate)) {
+                return std::unexpected("Fitting a radial guess requires a finite physical center.");
+            }
+        }
+        const long double stretch = static_cast<long double>(radius.value()) / seed.scales.radius.value();
+        // Apply each length factor to its physical quantity: s^5 alone can overflow while I*s^5 remains representable.
+        const double mass = static_cast<double>(seed.scales.mass.value() * stretch * stretch * stretch);
+        const double inertia = static_cast<double>(seed.moment_of_inertia * stretch * stretch * stretch * stretch * stretch);
+        const double rotation = static_cast<double>(seed.angular_velocity / stretch / stretch / stretch / stretch / stretch);
+        const double potential = static_cast<double>(center->potential.value() * stretch * stretch);
+        if (!std::isfinite(mass) || mass <= 0.0 || !std::isfinite(inertia) || inertia <= 0.0 ||
+            !std::isfinite(rotation) || (seed.angular_velocity != 0.0 && rotation == 0.0) || !std::isfinite(potential)) {
+            return std::unexpected("Fitting the radial guess exceeds representable mass, gravity or rotation scales. Use consistent physical units and a less extreme radius ratio.");
+        }
+        return FittedRadialSeed{.source = seed, .radius = radius, .mass = dimensions::MassValue{mass},
+            .moment_of_inertia = inertia, .angular_velocity = rotation};
+    }
+
+    /** @brief Evaluate unchanged material fields with the gravity generated by their stretched density. */
+    std::expected<utils::laneemden::PhysicalSample, std::string> evaluate_radial_seed(
+        const FittedRadialSeed& seed, dimensions::LengthValue radius) {
+        if (!std::isfinite(radius.value()) || radius.value() < 0.0 || !std::isfinite(seed.radius.value()) || seed.radius.value() <= 0.0 ||
+            !std::isfinite(seed.mass.value()) || seed.mass.value() <= 0.0 || !std::isfinite(seed.source.scales.radius.value()) || seed.source.scales.radius.value() <= 0.0 ||
+            !std::isfinite(seed.source.scales.gravitational_constant.value()) || seed.source.scales.gravitational_constant.value() <= 0.0) {
+            return std::unexpected("A fitted radial guess requires a finite nonnegative evaluation radius and finite positive source radius, fitted radius, mass and G.");
+        }
+        if (radius.value() >= seed.radius.value()) {
+            const long double gravitational_mass = static_cast<long double>(seed.source.scales.gravitational_constant.value()) * seed.mass.value();
+            const double potential = static_cast<double>(-gravitational_mass / radius.value());
+            const double gravity = static_cast<double>(gravitational_mass / radius.value() / radius.value());
+            if (!std::isfinite(potential) || !std::isfinite(gravity)) {
+                return std::unexpected("The fitted exterior gravity fields overflow in the supplied physical units.");
+            }
+            return utils::laneemden::PhysicalSample{
+                .density = dimensions::DensityValue{0.0}, .pressure = dimensions::PressureValue{0.0},
+                .enthalpy = dimensions::QuantityValue<dimensions::SpecificEnthalpy>{0.0},
+                .potential = dimensions::QuantityValue<dimensions::GravitationalPotential>{potential},
+                .enclosed_mass = seed.mass, .gravity = dimensions::QuantityValue<dimensions::Acceleration>{gravity}
+            };
+        }
+        const long double stretch = static_cast<long double>(seed.radius.value()) / seed.source.scales.radius.value();
+        const double source_radius = static_cast<double>((static_cast<long double>(radius.value()) / seed.radius.value()) * seed.source.scales.radius.value());
+        auto sample = evaluate_radial_seed(seed.source, dimensions::LengthValue{source_radius});
+        if (!sample) {
+            return sample;
+        }
+        const double mass = static_cast<double>(sample->enclosed_mass.value() * stretch * stretch * stretch);
+        const double potential = static_cast<double>(sample->potential.value() * stretch * stretch);
+        const double gravity = static_cast<double>(sample->gravity.value() * stretch);
+        if (!std::isfinite(mass) || !std::isfinite(potential) || !std::isfinite(gravity)) {
+            return std::unexpected("The fitted radial gravity fields overflow in the supplied physical units.");
+        }
+        sample->enclosed_mass = dimensions::MassValue{mass};
+        sample->potential = dimensions::QuantityValue<dimensions::GravitationalPotential>{potential};
+        sample->gravity = dimensions::QuantityValue<dimensions::Acceleration>{gravity};
+        return sample;
+    }
+}
diff --git a/src/lib/serif/solver/seed/polytrope.cpp b/src/lib/serif/solver/seed/polytrope.cpp
new file mode 100644
index 00000000..e2f35d7f
--- /dev/null
+++ b/src/lib/serif/solver/seed/polytrope.cpp
@@ -0,0 +1,60 @@
+#include "serif/solver/seed/polytrope.hpp"
+
+#include <cmath>
+#include <numbers>
+#include <utility>
+
+namespace serif::solver::seed {
+    std::expected<RadialSeed, std::string> make_polytropic_seed(
+        const eos::models::Polytrope& equation_of_state, dimensions::MassValue mass, PolytropeOptions options) {
+        const double n = equation_of_state.polytropic_index();
+        if (!std::isfinite(mass.value()) || mass.value() <= 0.0 || !std::isfinite(options.gravitational_constant.value()) || options.gravitational_constant.value() <= 0.0) {
+            return std::unexpected("A polytropic seed requires a positive finite target mass and gravitational constant.");
+        }
+        auto solution = utils::laneemden::solve(n, options.integration);
+        if (!solution) {
+            return std::unexpected(solution.error());
+        }
+        if (!solution->surface) {
+            return std::unexpected("A finite-radius stellar seed requires a Lane-Emden first zero; choose a finite-radius index and sufficient integration extent.");
+        }
+        if (!options.central_density) {
+            if (std::abs(n - 3.0) < 1.0e-12) {
+                return std::unexpected("An n=3 polytrope's mass is independent of central density. Specify central_density to choose its initial radius.");
+            }
+            const auto unit = utils::laneemden::scale_from_polytropic_constant(*solution, equation_of_state.polytropic_constant(),
+                dimensions::DensityValue{1.0}, options.gravitational_constant);
+            if (!unit) {
+                return std::unexpected(unit.error());
+            }
+            const double density = std::exp((2.0 * n / (3.0 - n)) * std::log(mass.value() / unit->mass.value()));
+            options.central_density = dimensions::DensityValue{density};
+        }
+        auto scales = utils::laneemden::scale_from_polytropic_constant(*solution, equation_of_state.polytropic_constant(),
+            *options.central_density, options.gravitational_constant);
+        if (!scales) {
+            return std::unexpected(scales.error());
+        }
+        // An explicit density defines an initial guess, including at n=3. Its physical mass
+        // need not equal the model's target; that prescription belongs to the equilibrium equations.
+
+        // Simpson integration of the regular dimensionless profile estimates I for the initial rigid rotation.
+        constexpr int intervals = 2048;
+        const double step = solution->surface->xi / intervals;
+        double integral = 0.0;
+        for (int index = 1; index < intervals; ++index) {
+            const double xi = index * step;
+            const auto sample = utils::laneemden::evaluate(*solution, xi);
+            if (!sample) {
+                return std::unexpected(sample.error());
+            }
+            const double weight = index % 2 == 0 ? 2.0 : 4.0;
+            integral += weight * std::pow(sample->theta, n) * std::pow(xi, 4);
+        }
+        const double inertia = (8.0 * std::numbers::pi / 3.0) * scales->central_density.value() * std::pow(scales->radial_scale, 5) * integral * step / 3.0;
+        if (!std::isfinite(inertia) || inertia <= 0.0) {
+            return std::unexpected("The spherical seed's moment of inertia is not finite and positive. Check the physical scales.");
+        }
+        return RadialSeed{.solution = std::move(*solution), .scales = *scales, .moment_of_inertia = inertia};
+    }
+}