Suppose we have $p_1,...,p_n$ and we have springs attached to them. We know that the lagrangian is $$ L = T - U = \sum_{i=1}^n \frac{1}{2} m_i \dot{p}_i^2 - \sum_{(i,j) \in E} \frac{1}{2} k_{ij} \left( \left\lVert p_i - p_j \right\rVert - r_{ij}\right)^2 $$
And the Euler lagrange equations are
$$ m_i \frac{d^2 p_i}{d t^2} = - \sum_{j \in \mathcal{N}(i)} k_{ij}\left( \left\lVert p_i - p_j \right\rVert - r_{ij}\right) \frac{p_i - p_j}{\left\lVert p_i - p_j \right\rVert} $$
Here $k_{ij}$ are elastic constants for the spring connecting vertex $i$ to vertex $j$, $E$ is the set of springs (modelled as pair of indices), $m_i$ is the i-th particle mass.
What I was wondering now suppose we have $p_k = p_k(t)$ for $k \in C \subset \left\{1,...,n \right\}$, so some of these particles have known trajectory (i.e. I am constraining them).
How do the Euler lagrange equation change? I've been trying to understand the how the constrained variations would apply in this case and I've been watching this tutorial but I struggle to work out the equation.
My only guess was to model everything as
$$ \begin{array}{ll} \min L = \min T - U & \\ s.t. \left\lVert p_k - p_k(t) \right\rVert^2 = 0 & \text{for $k \in C$} \end{array} $$
where $p_k(t)$ are given functions. I know in standard optimization we can use lagrange multipliers and it seems like according to the tutorial I shared we can extend the use to calculus of variations, I don't know why however I get confused when it comes to plug everything together.
Can anyone please clarify how to correctly include constraints in my case?
Update : This is what I have so far.. I am referring to Calculus of variations of Gelfand (Specifically theorem 2 of chapter 2, in section 12.2)
I rewrite my constraints as equivalent
$$ G(p_{k_1},\ldots,p_{k_N}) = \frac{1}{2} \sum_{k \in C} \left\lVert p_k - \hat{p}_k \right\rVert^2 $$
Where $\hat{p}$ are given functions. Using the theory of lagrange multipliers in this case I believe yields
$$ \frac{\partial L}{\partial q_i} + \lambda \frac{\partial G}{\partial q_i} - \frac{d}{dt} \frac{\partial L}{\partial \dot{q}_i} = 0 $$
Where given the type of constraint $\lambda$ is in general a function of $t$.
I have two cases to consider when deriving the system of differential equations if $i \notin C$ then we have $$ m_i \frac{d^2 p_i}{d t^2} = - \sum_{j \in \mathcal{N}(i)} k_{ij}\left( \left\lVert p_i - p_j \right\rVert - r_{ij}\right) \frac{p_i - p_j}{\left\lVert p_i - p_j \right\rVert} \;\;\;\;\; (1) $$
If $i \in C$ instead we have
$$ m_i \frac{d^2 p_i}{d t^2} = - \sum_{j \in \mathcal{N}(i)} k_{ij}\left( \left\lVert p_i - p_j \right\rVert - r_{ij}\right) \frac{p_i - p_j}{\left\lVert p_i - p_j \right\rVert} - \lambda (p_i - \hat{p_i}) $$
Because of my constraint $G$ however we have $p_i = \hat{p}_i$ so in truth we have
$$ \frac{d^2 \hat{p}_i}{d t^2} = - \sum_{j \in \mathcal{N}(i)} k_{ij}\left( \left\lVert \hat{p}_i - p_j \right\rVert - r_{ij}\right) \frac{\hat{p}_i - p_j}{\left\lVert \hat{p}_i - p_j \right\rVert} \;\;\;\;\; (2) $$
Now despite (2) and (1) might look the same what we have is a system of equation on $\left\{q_i\right\}$ however part of them are differential equations and part of them are algebraic equations (because if $i \in C$ then $\hat{p}_i$ is known). What I am looking for is for some numerical method that would enable me to solve the system of equations (differential + algebraic). I know methods for solving separately differential equations and algebraic equations, but I don't know if I wanted to combine them what I need to do.
Another alternative could be simplifying even further the system but then I don't really know how to do it.
Here I am now stuck because I don't know how to deal with this kind of systems. If anyone has a suggestion it would great.
Considering non trivial constraints of type $\|p_i-\hat{p}_i\|=c_i\gt 0$, which are holonomomic, those constraints can be integrated into the Lagrangian as follows:
$$ L = T-U + \sum_{j\in C}\lambda_j(\|p_j-\hat{p}_j\|-c_j) $$
so the movement equations are obtained from
$$ \cases{ m_i \frac{d^2 p_i}{d t^2} + \sum_{j \in \mathcal{N}(i)} k_{ij}\left( \left\lVert p_i - p_j \right\rVert - r_{ij}\right) \frac{p_i - p_j}{\left\lVert p_i - p_j \right\rVert}+\lambda_i\frac{p_i - \hat{p}_i}{\left\lVert p_i - \hat{p}_i \right\rVert}=0\\ \|p_j-\hat{p}_j\|-c_j=0 } $$
so we have $\lambda_j,\ \ j\in C$ additional unknowns and $\|p_j-\hat{p}_j\|-c_j=0,\ \ j\in C$ additional equations.
NOTE
To avoid the solution of a system of differential-algebraic equations, we can consider instead $\frac{d^2}{dt^2}\left(\|p_j-\hat{p}_j\|-c_j\right)=0$ and then we can obtain a set of ordinary differential equations for $\frac{d^2p_i}{dt^2},\ \ i\in N$, which are $\lambda_j$ independent, being aware that the initial conditions for this ODE system should be consistent with the constraints.
Attached a MATHEMATICA script modeling the reduced case study with three points in the plane, with one of then constrained by $\|p_1-\hat{p}_1\|= c_1$ where $\hat{p}_1$ describes a centered circle with radius $r$ and angular velocity $\omega$
Follows a plot showing in red blue and green respectively, the springs $(p_1,p_2)$, $(p_2,p_3)$ and $(p_1,p_3)$.
Follows also another plot showing the preservation of distance $c_1$ between $p_1$ and $q_1$ representing $\hat{p}_1$.