This is my first question on stackoverflow. Sorry if it is a bit clumsy.
I am trying to solve a static equilibrium describing a beam in a layered soil. In terms of forces acting on the beam, there is some external force, and there is the reaction of the soil to that force. The action and reaction obviously have to balance out. The reaction is the sum of forces from the soil integrated over the different layers. The soil responds non-linearly as a function of displacement. Basically, this particular equilibrium is a scalar because it's just the sum of reaction = action.
On the other hand, the internal forces within the beam also have to balance with the local force acting on the beam element in question. Let's say we have 40 such elements. So the problem I face in a nutshell is the following:
I need to find the root of F1(x1,x2,...x40), with x1,x2...x40 describing the displacement of the pile in the soil. I pass and initial guess for the displacement, and Fsolve comes back with a vector length 40 which gives F1 = [0,0...0].. So this basically balances the internal forces with the local force on the beam (actually the curvature with moment).
The problem is I cannot seem to add an additional constraint which is the global equilibrium that I described above because that violates the condition in fsolve that the output vector length has to equal the input vector length. I cannot seem to add such constraint within Fsolve. I have also gone the way of optimize where you can add constraints, but that doesn't work because it optimizes a scalar, and does not work with vector roots.
Perhaps I just missed something stupid, but I've scoured stackoverflow in vain a few days. I can probably solve it with differential equation solvers but i don't want to go that way for various reasons if I can avoid it. Anything to give me a new angle would help! Brand new to python, so I much appreciate all the helpful solutions I've read so far even if it didn't quite get me there.
Usually, we have n equations and n variables, i.e. a square system. If the solution is not unique, it means that your equations are not independent. When adding additional constraints, you no longer can use a square system solver. However, it is quite possible to use a nonlinear optimization solver. E.g. we can write:
min 0
subject to
F(x) = 0 (n equations)
G(x) = 0 (additional constraints)
If you want a solution that is close, you can add slacks and minimize the sum of the squared errors:
min ||s|| + ||t||
F(x) = s
G(x) = t
I use this approach quite a lot.