homotopy method and warm start

[tvdbogert]

Huawei Wang has modified the double inverted pendulum notebook to do control identification on actual human data. He finds a different local optimum after each restart with random initial guess. The controller model at this time is simple and far from human, but even with a bad model, the parameter identification should find the unique global optimum.

We want to implement the homotopy method from the Vyasarayani paper. This requires:

1. An extra term in the implicit dynamics which tracks the data. It would be easy to do in the notebook, but we are wondering if this is general enough that it can be added as an option in the opty toolbox.
2. Doing a sequence of IPOPT runs while the lambda parameter increases. We have figured out how to warm start IPOPT (in Matlab) with lagrange multipliers from a previous runs and this works well. Ideally this would all be done inside opty.

@moorepants please comment on these feature requests and give some guidance on where and how to do this. Huawei is pretty comfortable with python, we just need some guidance.

[moorepants]

This all sounds great. I think adding the features so that it is general sounds fine. So, if Huawei is comfortable with the code he can just start adding things and then submit a pull request and we can discuss the changes. Of if he wants to call me and describe what specifically he wants to do, we can make a plan before he starts writing code. Seeing some examples of what you have in mind on a simple problem could help me think about how to generalize it. Glad you guys are trying this out!

[HuaweiWang]

@moorepants Thanks a lot for your help.

The simple example of what we want to do is as follow:

Problem 1: adding extra term
In the original problem, the implicit dynamics is f(x, x_dot, p, r) = 0. Well, by applying homotopy method, we need to reverse implicit dynamics to (1- lamda)f(x, x_dot, p, r) + Lamda*(x_dot – (x_data – x)) = 0, in which Lamda will decrease from 1 to 0.

We want add the extra part, Lamda*(x_dot – (x_data – x)), into opty codes to have the option that if we want use homotopy method, then the implicit dynamics will have the extra part; otherwise the implicit dynamics will still be f(x, x_dot, p, r) = 0.

Problem 2: warm start
We want to active warm start function in IPOPT to achieve the homotopy method. But we don't know the setting commands in python IPOPT right now.

As I haven't read the code in opty before, I think it will be better for us to discuss after I am familiar with the opty. I will try to add things first and look for your help if I got stuck!

[moorepants]

In the original problem, the implicit dynamics is f(x, x_dot, p, r) = 0. Well, by applying homotopy method, we need to reverse implicit dynamics to (1- lamda)f(x, x_dot, p, r) + Lamda*(x_dot – (x_data – x)) = 0, in which Lamda will decrease from 1 to 0.

Ok, that seems easy enough. Let me know if you need advice on how this should be added.

We want to active warm start function in IPOPT to achieve the homotopy method. But we don't know the setting commands in python IPOPT right now.

The settings in Python are identical to those you give in C/C++.

[HuaweiWang]

OK, I will let you know if I need help. thanks.

[tvdbogert]

Huawei, I looked at the code you just posted. Some quick comments:

• typing error: Lamda and Lamada should both be Lambda
• I noticed that you only added the "data tracking" dynamics term to equations 3 and 4. Vyasarayani et al. have this term in all equations. I like your version because equations 1 and 2 are just integrators and do not respond to the surface accelerations. However, that makes the tracking dynamics a second order system without damping. You will probably need to add damping or put the tracking term also in equations 1 and 2.
• The "data tracking" term also needs a gain K. Possibly K=1 is OK. In equations 3 and 4, K will have units of 1/s^2 and is roughly the inverse of the square of the time constant of the tracking dynamics (when lambda=1).

[HuaweiWang]

Hi Jason,

Thanks for the comments.

• Response to comment 1 : I will change them, thanks.

• Response to comment 2: the "data tracking" dynamics term we used for
  equations 3 and 4 is as follow:

  (1-Lambda)* f3 + Lambda_(Omiga_a_dot - (Omiga_a_data - Omiga_a)) = 0
  (1-Lambda)_f4 + Lambda*(Omiga_h_dot - (Omiga_h_data - Omiga_h)) = 0

when Lambda = 1, we got Omiga_a = Omiga_a_data/(s+1) and Omiga_h =
Omiga_h_data/(s+1).

if we do not apply "data tracking" dynamics term on equations 1 and 2, then
Theta_a = Omiga_a_data/(s_(s+1)) and Theta_h = Omiga_h_data/(s_(s+1)).
if we apply "data tracking" dynamics term on equations 1 and 2, then
Theta_a = Theta_a_data/(s+1) and Theta_h = Theta_h_data/(s*(s+1)).
As in our experiment data, Theta_data and Omiga_data has a integral
relationship, I think in both ways, we will get the same result.

• Response to comment 3: you are right, we need a K to make sure the
  "data tracking "dynamic term faster than actual model (experiment data
  dynamics).

Bests,

Huawei

On Thu, May 19, 2016 at 12:57 PM, Ton van den Bogert <
notifications@github.com> wrote:

Huawei, I looked at the code you just posted. Some quick comments:

• typing error: Lamda and Lamada should both be Lambda
• I noticed that you only added the "data tracking" dynamics term to
  equations 3 and 4. Vyasarayani et al. have this term in all equations. I
  like your version because equations 1 and 2 are just integrators and do not
  respond to the surface accelerations. However, that makes the tracking
  dynamics a second order system without damping. You will probably need to
  add damping or put the tracking term also in equations 1 and 2.
• The "data tracking" term also needs a gain K. Possibly K=1 is OK. In
  equations 3 and 4, K will have units of 1/s^2 and is roughly the inverse of
  the square of the time constant of the tracking dynamics (when lambda=1).

—
You are receiving this because you are subscribed to this thread.
Reply to this email directly or view it on GitHub
#39 (comment)

[moorepants]

@HuaweiWang I think we should do a google hangout for about a 30 min to an hour to discuss exactly what you want to do and I can give you some guidance on how to implement it in the code. I think it will be best to make a basic plan first.

[HuaweiWang]

Sure, what time do you available?

I am available now to 5:00pm

On Thu, May 19, 2016 at 2:34 PM, Jason K. Moore notifications@github.com
wrote:

@HuaweiWang https://github.com/HuaweiWang I think we should do a google
hangout for about a 30 min to an hour to discuss exactly what you want to
do and I can give you some guidance on how to implement it in the code. I
think it will be best to make a basic plan first.

—
You are receiving this because you were mentioned.
Reply to this email directly or view it on GitHub
#39 (comment)

[moorepants]

I can speak right now for about 30 minutes. My gmail account is moorepants at gmail.

[HuaweiWang]

Hi Jason,

Do you received my invitation?

On Thu, May 19, 2016 at 3:09 PM, Jason K. Moore notifications@github.com
wrote:

I can speak right now for about 30 minutes. My gmail account is moorepants
at gmail.

—
You are receiving this because you were mentioned.
Reply to this email directly or view it on GitHub
#39 (comment)