1
SFML projects / Re: Drawing a 3d world using SFML and a perspective tranformation.
« on: February 10, 2024, 08:34:43 pm »
Hi all!
I'm adding this post for those who are curious about how I handled the motion of the pendulum in the 1st video. It goes to the limits of my abilities in several ways.
I had an opportunity while working 3rd shift for a couple of years to finally
take on that long overdue review of classical mechanics. I had (typically) 6 hours per shift of free time, and it was nice being paid $22/hour to do it.
I obtained and studied from: Introduction to Classical Mechanics by David Morin. 1st edition (in case there's a 2nd by now).
My educational background: BS in physics long ago.
I'm using skills which would have faded to nothing by now if not for my usage in this hobby.
My goal was to make it through chapter 6 on the Lagrangian method. It's an alternative formulation (to vector based Newtonian mechanics) based on scalar quantities (just a number) such as kinetic and potential energy and is well suited to problems where motion is limited by a constraint, such as being limited to motion on the surface of a sphere in 3d space like this pendulum.
I highly recommend the textbook. It is as hardcore as anyone could want and it's loaded with fully worked solutions for a portion of the very full range of problems in each chapter. It's also loaded with clever little physics limerics. I'll post 1 or 2 here later.
This text has a place of honor in my work here. It can be seen floating in midair in the room in the 2nd video.
Please find attached type persPendulum. It's virtual update() has a bit more than usual going on, which gets right away to an issue with the enormous time step between frames, and to why type vec3d became involved.
My simplistic numerical integration scheme is this:
step 1: obtain the acceleration to be applied this frame.
In the case of the shots the acceleration is constant = (0,g,0).
In the case of the pendulum it's far from that.
step 2:
update the position:
Then update the velocity
The analysis will be in spherical coordinates.
So far these methods are persPendulum methods only, but I think they also
should be static methods available in vec3d?
I'll also attach my solution on paper soon, but I want it to be much cleaner than the working version I have.
The identifiers for the two angles theta (polar) and phi (azimuthal) and their derivatives are:
// polar angle and 1st and 2nd time derivatives
double th, th_1, th_2;
// azimuthal angle and 1st and 2nd time derivatives
double phi, phi_1, phi_2;
The 2 equations of motion obtained by the Lagrangian method are:
An additional condition may also be applied, perhaps to stabilize the integration?
Because the Lagrangian = kinetic - potential energy doesn't depend explicitly on the angle phi, there is a constant of the motion. Here it's the z component of the angular momentum so:
So the units of measure ought to balance above.
Integration at the frame rate ( dt = 20ms ) does not give a stable result.
It became stable at 10 updates per frame, so I tripled that and called it good.
The pendulum in the video is being updated 30 times per frame.
Once all of these micro updates are done I assign bobs position and velocity.
I'm still struggling with the float/double interface here.
I'll leave off here for now. But I must deal with the sudden change in th_1 and phi_1 which occur when poor bob gets shot.
I'm back to wrap up.
As for handling bobs hits:
I'm cheating a bit. I call hitFree() and do not account for the constraint in the collision,
so the momentum transferred isn't accurate, but it seems to do fine for effect. This serves the present purpose well enough for me. Next problem please.
Update the values for th_1 and phi_1
bobs velocity in polar coordinates is:
Taking the dot product (of the return value above) with thHat (unit vector in theta direction) gives
Vel.dot( thHat ) = L*th_1
solving for th_1 permits the update
Just call
Why dt? It won't be moved here, just in update(). The collision test is ray based so dt is used to obtain the position in the previous frame = 2 points for a ray.
Finally, I'll admit that the integration details are well beyond my knowledge.
I think I should be using th_1 = Lz/( sin(th)*sin(th) ); instead of the 2nd acceleration equation, which actually follows from above. I'm presently experimenting with using
However, I'm once again getting runaway cases...
It's a fight to maintain stability and I don't know enough right here.
Perhaps someone who knows this part better may advise? Please?
I've also now attached my "on paper" solution.
And that's all I've got so I'll call it a wrap.
Any questions and/or suggestions are welcome.
I'm adding this post for those who are curious about how I handled the motion of the pendulum in the 1st video. It goes to the limits of my abilities in several ways.
I had an opportunity while working 3rd shift for a couple of years to finally
take on that long overdue review of classical mechanics. I had (typically) 6 hours per shift of free time, and it was nice being paid $22/hour to do it.
I obtained and studied from: Introduction to Classical Mechanics by David Morin. 1st edition (in case there's a 2nd by now).
My educational background: BS in physics long ago.
I'm using skills which would have faded to nothing by now if not for my usage in this hobby.
My goal was to make it through chapter 6 on the Lagrangian method. It's an alternative formulation (to vector based Newtonian mechanics) based on scalar quantities (just a number) such as kinetic and potential energy and is well suited to problems where motion is limited by a constraint, such as being limited to motion on the surface of a sphere in 3d space like this pendulum.
I highly recommend the textbook. It is as hardcore as anyone could want and it's loaded with fully worked solutions for a portion of the very full range of problems in each chapter. It's also loaded with clever little physics limerics. I'll post 1 or 2 here later.
This text has a place of honor in my work here. It can be seen floating in midair in the room in the 2nd video.
Please find attached type persPendulum. It's virtual update() has a bit more than usual going on, which gets right away to an issue with the enormous time step between frames, and to why type vec3d became involved.
My simplistic numerical integration scheme is this:
step 1: obtain the acceleration to be applied this frame.
In the case of the shots the acceleration is constant = (0,g,0).
In the case of the pendulum it's far from that.
step 2:
update the position:
position += velocity*dt + 0.5f*acceleration*dt*dt;
This is from Taylors series to 2 terms. Then update the velocity
velocity += acceleration*dt;
The analysis will be in spherical coordinates.
So far these methods are persPendulum methods only, but I think they also
should be static methods available in vec3d?
// base vectors in spherical coordinates
vec3d get_rHat()const{ return sin(th)*( cos(phi)*Xu + sin(phi)*Yu ) - cos(th)*Zu; }
vec3d get_thHat()const{ return cos(th)*( cos(phi)*Xu + sin(phi)*Yu ) + sin(th)*Zu; }
vec3d get_phiHat()const{ return cos(phi)*Yu - sin(phi)*Xu; }
Where the vectors Xu, Yu, Zu are the cartesian basis.vec3d get_rHat()const{ return sin(th)*( cos(phi)*Xu + sin(phi)*Yu ) - cos(th)*Zu; }
vec3d get_thHat()const{ return cos(th)*( cos(phi)*Xu + sin(phi)*Yu ) + sin(th)*Zu; }
vec3d get_phiHat()const{ return cos(phi)*Yu - sin(phi)*Xu; }
I'll also attach my solution on paper soon, but I want it to be much cleaner than the working version I have.
The identifiers for the two angles theta (polar) and phi (azimuthal) and their derivatives are:
// polar angle and 1st and 2nd time derivatives
double th, th_1, th_2;
// azimuthal angle and 1st and 2nd time derivatives
double phi, phi_1, phi_2;
The 2 equations of motion obtained by the Lagrangian method are:
// step 1: find acceleration
th_2 = phi_1*phi_1*sinTh*cosTh - GravY*sinTh/L;// for theta
phi_2 = -2.0*th_1*phi_1*cosTh/sinTh;// for phi
Here, sinTh, cosTh are saved values of sin(th), cos(th), GravY = gravity.y and L = length of support. Note that the mass isn't a factor.th_2 = phi_1*phi_1*sinTh*cosTh - GravY*sinTh/L;// for theta
phi_2 = -2.0*th_1*phi_1*cosTh/sinTh;// for phi
An additional condition may also be applied, perhaps to stabilize the integration?
Because the Lagrangian = kinetic - potential energy doesn't depend explicitly on the angle phi, there is a constant of the motion. Here it's the z component of the angular momentum so:
phi_1 = Lz/( sin(th)*sin(th) );// based on Lz conservation
where the value of Lz is determined by initial conditions. The quantity here has been scaled from the actual Lz: Scaled this is actually Lz /= m*L*L.So the units of measure ought to balance above.
Integration at the frame rate ( dt = 20ms ) does not give a stable result.
It became stable at 10 updates per frame, so I tripled that and called it good.
The pendulum in the video is being updated 30 times per frame.
Once all of these micro updates are done I assign bobs position and velocity.
I'm still struggling with the float/double interface here.
// bobs new velocity
vec3d thHat = cosTh*( cosPhi*Xu + sinPhi*Yu ) + sinTh*Zu;// current thHat
vec3d phiHat = cosPhi*Yu - sinPhi*Xu;// current phiHat
Vel = L*( th_1*thHat + sinTh*phi_1*phiHat );
bob.vel = Vel.to_vec3f();// bob.vel update
vec3d rHat = sinTh*( cosPhi*Xu + sinPhi*Yu ) - cosTh*Zu;// updated rHat
// and position
vec3d Pos = L*rHat;
Pos.x += static_cast<double>(topQ.pos.x);
Pos.y += static_cast<double>(topQ.pos.y);
Pos.z += static_cast<double>(topQ.pos.z);
bob.setPosition( Pos.to_vec3f() );// bob update
vec3d thHat = cosTh*( cosPhi*Xu + sinPhi*Yu ) + sinTh*Zu;// current thHat
vec3d phiHat = cosPhi*Yu - sinPhi*Xu;// current phiHat
Vel = L*( th_1*thHat + sinTh*phi_1*phiHat );
bob.vel = Vel.to_vec3f();// bob.vel update
vec3d rHat = sinTh*( cosPhi*Xu + sinPhi*Yu ) - cosTh*Zu;// updated rHat
// and position
vec3d Pos = L*rHat;
Pos.x += static_cast<double>(topQ.pos.x);
Pos.y += static_cast<double>(topQ.pos.y);
Pos.z += static_cast<double>(topQ.pos.z);
bob.setPosition( Pos.to_vec3f() );// bob update
I'll leave off here for now. But I must deal with the sudden change in th_1 and phi_1 which occur when poor bob gets shot.
I'm back to wrap up.
As for handling bobs hits:
I'm cheating a bit. I call hitFree() and do not account for the constraint in the collision,
so the momentum transferred isn't accurate, but it seems to do fine for effect. This serves the present purpose well enough for me. Next problem please.
Update the values for th_1 and phi_1
bobs velocity in polar coordinates is:
vec3d persPendulum::get_vel()const
{ return vec3d( L*( th_1*get_thHat() + sin(th)*phi_1*get_phiHat() ) ); }
Let vec3d Vel = bob.vel following the collision, then{ return vec3d( L*( th_1*get_thHat() + sin(th)*phi_1*get_phiHat() ) ); }
Taking the dot product (of the return value above) with thHat (unit vector in theta direction) gives
Vel.dot( thHat ) = L*th_1
solving for th_1 permits the update
th_1 = Vel.dot(thHat)/L;
Take the dot product with phiHat instead to obtainphi_1 = Vel.dot( phiHat )/( L*sin(th) );
following this I update Vel and then bob.velVel = get_vel();
bob.vel = Vel.to_vec3f();
bob.vel = Vel.to_vec3f();
Just call
bool hit( persBall& PB, float cR, float dt );
where Cr = coefficient of restitution allows partially elastic collisions.Why dt? It won't be moved here, just in update(). The collision test is ray based so dt is used to obtain the position in the previous frame = 2 points for a ray.
Finally, I'll admit that the integration details are well beyond my knowledge.
I think I should be using th_1 = Lz/( sin(th)*sin(th) ); instead of the 2nd acceleration equation, which actually follows from above. I'm presently experimenting with using
phi_1 = Lz/( sinTh*sinTh );
phi += phi_1*dT;
in place ofphi += phi_1*dT;
phi_2 = -2.0*th_1*phi_1*cosTh/sinTh;
phi += phi_1*dT + 0.5*phi_2*dT*dT;
phi_1 += phi_2*dT;
and it's giving a stable result for lower # of micro updates.phi += phi_1*dT + 0.5*phi_2*dT*dT;
phi_1 += phi_2*dT;
However, I'm once again getting runaway cases...
It's a fight to maintain stability and I don't know enough right here.
Perhaps someone who knows this part better may advise? Please?
I've also now attached my "on paper" solution.
And that's all I've got so I'll call it a wrap.
Any questions and/or suggestions are welcome.