FAQ  •  Login

Matlab Assignment 9

<<

F11Carrie

Newbie
Newbie

Posts: 23

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Fri Dec 02, 2011 1:10 pm

Matlab Assignment 9

I'm a little lost on the Newmark integration method.

1. From what I understand the first step is to use Newton's Equations of Motion just the way we have been to find qdotdot and lambda at time t=0, so that we know everything about the system at that time.

2. Then we start incrementing time. At each time step our first guess for qdotdot and lambda is the value from the previous time step and we use these values to guess updated values of q and qdot based on Newmark's equations. Based on this guess of q, we evaluate J-tilde - a value we'll use/save for every iteration in this time step.

3.From here we can use the quasi-Newton method to get a more accurate value for qdotdot and lambda using
Psi=[M*qdotdot+[Phiq]T(current_guess_qn+1)-QA(qdot,q,t);
1/(Beta*h2)*Phi(current_guess_qn+1,t)]

4. Then we update our guesses for q and qdot, check the size of the correction factor and decide whether to go through the quais-Newton again to get better values for qdotdot and lambda (back to 3 if correction is still larger than specified limit).

5. Once we're satisfied with qdotdot and lambda, we save those values along with q and qdot computed from those values and then move on to the next time step and repeat the process starting at 2.
<<

Dan Negrut

Global Moderator
Global Moderator

Posts: 833

Joined: Wed Sep 03, 2008 12:24 pm

Unread post Fri Dec 02, 2011 6:54 pm

Re: Matlab Assignment 9

F11Carrie wrote:I'm a little lost on the Newmark integration method.

1. From what I understand the first step is to use Newton's Equations of Motion just the way we have been to find qdotdot and lambda at time t=0, so that we know everything about the system at that time.

2. Then we start incrementing time. At each time step our first guess for qdotdot and lambda is the value from the previous time step and we use these values to guess updated values of q and qdot based on Newmark's equations. Based on this guess of q, we evaluate J-tilde - a value we'll use/save for every iteration in this time step.

3.From here we can use the quasi-Newton method to get a more accurate value for qdotdot and lambda using
Psi=[M*qdotdot+[Phiq]T(current_guess_qn+1)-QA(qdot,q,t);
1/(Beta*h2)*Phi(current_guess_qn+1,t)]

4. Then we update our guesses for q and qdot, check the size of the correction factor and decide whether to go through the quasi-Newton again to get better values for qdotdot and lambda (back to 3 if correction is still larger than specified limit).

5. Once we're satisfied with qdotdot and lambda, we save those values along with q and qdot computed from those values and then move on to the next time step and repeat the process starting at 2.


Carrie - Actually, this is it.
One more thing: don't look at the value of \delta lambda to understand whether you need to keep iterating or not. Rather, take a look at the infinity norm of \delta qdotdot. If it's less then 1E-3, you can bail out of the quasi-Newton.
I hope this helps.
Dan
<<

F11Jonathan

Newbie
Newbie

Posts: 29

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Sat Dec 03, 2011 3:22 pm

Re: Matlab Assignment 9

My solution for k=2 is unstable. The acceleration increases to infinity. Do I need an initial guess for acceleration and lambda? I am using zeros for now.

The order of operations in my loop are:
while norm>1*10^3
1.recalculate q(n+1) and q_d(n+1) using updated q_dd(n+1)
2.recalculate phi and the jacobian using the newest values for q(n+1) q_d(n+1) and q_dd(n+1)
3.J tilde calculated only in the first iteration (or in every iteration, I've tried both)
4.Calculate the residual of q_dd and lambda
5.Calculate the correction
6.Update q_dd(n+1) and lambda(n+1) by subtracting the correction
7.Calculate the norm of the q_dd correction
end

with lambda(1)=[0 0 0 0]'
and q_dd(1)=[0 0 0 0]'
after the first iteration, my residual is
[0 19.62 -2.5 -5 19.62 0 0 0 0 0]'
which is -Qa and zeros for the scaled phi
after the first iteration, my correction is
[0 39.24 -1.6675 -10 39.24 0 0 22.12 5 -2.5]'

For the Newmark method my beta is 0.3025 and gamma is 0.6
<<

Dan Negrut

Global Moderator
Global Moderator

Posts: 833

Joined: Wed Sep 03, 2008 12:24 pm

Unread post Sat Dec 03, 2011 4:10 pm

Re: Matlab Assignment 9

F11Jonathan wrote:My solution for k=2 is unstable. The acceleration increases to infinity. Do I need an initial guess for acceleration and lambda? I am using zeros for now.

The order of operations in my loop are:
while norm>1*10^3
1.recalculate q(n+1) and q_d(n+1) using updated q_dd(n+1)
2.recalculate phi and the jacobian using the newest values for q(n+1) q_d(n+1) and q_dd(n+1)
3.J tilde calculated only in the first iteration (or in every iteration, I've tried both)
4.Calculate the residual of q_dd and lambda
5.Calculate the correction
6.Update q_dd(n+1) and lambda(n+1) by subtracting the correction
7.Calculate the norm of the q_dd correction
end

with lambda(1)=[0 0 0 0]'
and q_dd(1)=[0 0 0 0]'
after the first iteration, my residual is
[0 19.62 -2.5 -5 19.62 0 0 0 0 0]'
which is -Qa and zeros for the scaled phi
after the first iteration, my correction is
[0 39.24 -1.6675 -10 39.24 0 0 22.12 5 -2.5]'

For the Newmark method my beta is 0.3025 and gamma is 0.6


Use the values from the previous time step to seed qdotdot and lambda. At time t=0 solve the EOM in conjunction with the Acc Kin Constr Eq to get qdotdot_0 and \lambda_0.
Also, in "4.Calculate the residual of q_dd and lambda", the residual is the value of Psi, not clear exactly what you mean by "residual of q_dd and lambda".
Make sure your step-size h is small enough.
How many constraints and how many generalized coordinates do you have? Do you pose the problem right?
Dan
<<

F11Jonathan

Newbie
Newbie

Posts: 29

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Sat Dec 03, 2011 4:33 pm

Re: Matlab Assignment 9

It's working now

At time t=0 solve the EOM in conjunction with the Acc Kin Constr Eq to get qdotdot_0 and \lambda_0

I was just using zeros for the initial guess
<<

F11Wendy

Newbie
Newbie

Posts: 26

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Sun Dec 04, 2011 9:15 am

Re: Matlab Assignment 9

Question 1
#1 Plots.jpg
#1 Plots.jpg (38.08 KiB) Viewed 8352 times
<<

F11Christopher

Newbie
Newbie

Posts: 22

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Mon Dec 05, 2011 12:34 am

Re: Matlab Assignment 9

I'm stuck on how to use newton-raphson in conjugation with the backwards Euler method in problem 2. Once we run our initial guess through the equation, how do we calculate the correction factor? I thought we just divide our answer by the derivative of the equation, but this is giving me messy results. Any pointers?
<<

Dan Negrut

Global Moderator
Global Moderator

Posts: 833

Joined: Wed Sep 03, 2008 12:24 pm

Unread post Mon Dec 05, 2011 8:17 am

Re: Matlab Assignment 9

F11Christopher wrote:I'm stuck on how to use newton-raphson in conjugation with the backwards Euler method in problem 2. Once we run our initial guess through the equation, how do we calculate the correction factor? I thought we just divide our answer by the derivative of the equation, but this is giving me messy results. Any pointers?


Let's say that you are at t0 and want to find y1 at t1. Then you discretize your ODE and end up with an equation that looks like \Phi(y1)=y1 - h*f(t1,y1)=0.
Now you use Newton Raphson to find y1.
a) Take an initial guess (value of y at previous time step should be good) for y1 and set t = t+h.
b) Compute the residual corresponding to the current value of y1
c) Compute the Jacobian (sensitivity) \Psi_y1 = 1- h*partial f/partial y
d) Compute the correction \delta y = residual/Psi_y1
e) Apply the correction \delta y: y_improved = y - \delta y.
f) If correction small enough move on to the next time step (go to a) above). If not, go to b) above.

I think that the most important step is to understand the form of the nonlinear equation \Phi(y1)=y1 - h*f(t1,y1)=0 and where it comes from. After that you simply apply Newton-Raphson.
I hope this helps.
Dan
<<

F11Jonathan

Newbie
Newbie

Posts: 29

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Mon Dec 05, 2011 1:24 pm

Re: Matlab Assignment 9

Here is what I have so far.
Problem 1 Solution
Problem1.jpg
Problem1.jpg (15.66 KiB) Viewed 8309 times

Problem 1 between solution and ODE 45 solution
Problem1Bonus.jpg
Problem1Bonus.jpg (15.17 KiB) Viewed 8309 times

Problem 2 Forward Euler and Backward Euler Solutions
Problem2-1.jpg
Problem2-1.jpg (14.54 KiB) Viewed 8309 times
<<

F11Jonathan

Newbie
Newbie

Posts: 29

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Mon Dec 05, 2011 1:27 pm

Re: Matlab Assignment 9

Problem 2 Difference Between Forward and Backward Euler Solutions
Problem2-2.jpg
Problem2-2.jpg (15.03 KiB) Viewed 8309 times

Problem 3 Order Analysis
Problem3.jpg
Problem3.jpg (13.9 KiB) Viewed 8309 times

Problem 4 BDF Solution
Problem4.jpg
Problem4.jpg (13.79 KiB) Viewed 8309 times
<<

F11Carrie

Newbie
Newbie

Posts: 23

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Mon Dec 05, 2011 3:35 pm

Re: Matlab Assignment 9

Problems 1 and 2.jpg
My solutions for the Forward and Backward Euler Methods
Problems 1 and 2.jpg (93.16 KiB) Viewed 8303 times

Problem_3.jpg
Problem 3 - Order analysis of Forward Euler vs. Runge-Kutta
Problem_3.jpg (76.87 KiB) Viewed 8303 times

Problem_4.jpg
Problem 4 - 2nd Order BDF
Problem_4.jpg (21.21 KiB) Viewed 8303 times
<<

F11Carrie

Newbie
Newbie

Posts: 23

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Mon Dec 05, 2011 3:37 pm

Re: Matlab Assignment 9

Problem_5.jpg
Problem 5 - Newmark integration results
Problem_5.jpg (132.8 KiB) Viewed 8302 times


Note that torque produced by the revolute joint at the location of the joint is zero because the reaction forces pass through the point where torque is calculated. This value would be non-zero if the torque was calculated at the centroid.
<<

F11Jonathan

Newbie
Newbie

Posts: 29

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Mon Dec 05, 2011 4:48 pm

Re: Matlab Assignment 9

Here is Problem 5 Force and Torque
Problem5-F.jpg
Problem5-F.jpg (19.52 KiB) Viewed 8297 times

Problem5-T.jpg
Problem5-T.jpg (13.29 KiB) Viewed 8297 times

It looks like they agree with Carrie's
<<

F11Wendy

Newbie
Newbie

Posts: 26

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Mon Dec 05, 2011 9:31 pm

Re: Matlab Assignment 9

#1 Plots.jpg
#1 Plots.jpg (38.08 KiB) Viewed 8286 times

#2 Plots.jpg
#2 Plots.jpg (37.39 KiB) Viewed 8286 times

Prob 3.jpg
Prob 3.jpg (36.74 KiB) Viewed 8286 times
<<

F11Wendy

Newbie
Newbie

Posts: 26

Joined: Mon Sep 12, 2011 7:58 pm

Unread post Mon Dec 05, 2011 9:32 pm

Re: Matlab Assignment 9

#4 Plot.jpg
#4 Plot.jpg (18.83 KiB) Viewed 8286 times

Prob 5.jpg
Prob 5.jpg (50.68 KiB) Viewed 8286 times
Next

Return to ME451 Fall 2011: Kinematics and Dynamics of Machine Systems

Who is online

Users browsing this forum: No registered users and 1 guest

Powered by phpBB © 2000, 2002, 2005, 2007 phpBB Group.
Designed by ST Software.