Function will not converge

[remslie]

Hi,

I am calculating dynamic pressures on a moving gate and cannot get the integral for vertical force Fv(t) to converge. I have tried adjusting the convergence tolerance but this does not seem to make any difference. Any clues would be appreciated.

Cheers

Ross

[AlanStevens]

See if the attached helps.

Alan

[remslie]

Alan,

Yes indeed your modification did help. If I now take this function and wish to use it in an Odesolve block as the final step I hit the wall with another error. Any clues to de-bug this would be appreciated.

The problem I am trying to solve can be unstable with divergent response, however, it should be stable and converge when the damping (zeta) is a modest value (i have used 0.05 or 5%).

Thanks again

ross

[AlanStevens]

Ross

The attached works - the ODE solver worked ok, it was just the plotting that confused it (because of a previously defined set of t's).

However, I notice that the solve block contains M(t), which is a function of Fv(t), which is in turn a function of pp(t,y), which is a function of dtheta/dt - in your case a previously set function of dtheta/dt, not the one you are currently solving for within the solve block. Is this what you intend?

Alan

Edit. Also, it is a function of theta, through F(t,s). Again this isn't the theta that the ODE block is solving for. Is this right?

[Werner_E]

Edit. Also, it is a function of theta, through F(t,s). Again this isn't the theta that the ODE block is solving for. Is this right?

Guess no and it makes quite a difference. See below.

BTW, I wondered why everytime I ran the sheet I had an Excel sheet ("height 2.xls") in the same directory the file would be in. At last I found that the evaluation of Z= is the cause - in its properties the creation of that file is activated. It was first introduced in Ross' version 5 sheet. Nasty, that hidden way to write to my disk. Its cleared in the attached sheet.

[AlanStevens]

Although this gets in the thetadot term it doesn't account for the theta term under the integral in F(t,s), which is part of the definition of pp.

Hadn't noticed the spreadsheet! Good catch!

Alan

[Werner_E]

AlanStevens wrote:

Although this gets in the thetadot term it doesn't account for the theta term under the integral in F(t,s), which is part of the definition of pp.

Oh,yes. I missed that. Was a quick hack - too quick! So the the sheet would need some more modifications as I think odesolvbe will not work with a function to solve for as a parameter of another function involved.

Hadn't noticed the spreadsheet! Good catch!

To be honest I wasn't really fully aware that its so easy to export data using a simple evaluation 😉

[AlanStevens]

Werner Exinger wrote:

So the the sheet would need some more modifications as I think odesolvbe will not work with a function to solve for as a parameter of another function involved.

I think one would have to obtain dF(t,s)/dt and solve for F simultaneously with theta etc. However, it's possible that Ross has just used the same symbol (theta) for two entirely different parameters, so it's not worth thinking about in more detail until we hear from him!

Alan

[remslie]

Alan,Werner,

Again thanks for your interest.

Theta(t) is indeed consistent in all of the integral equations and equation of motion in the solve block.We could delete the assignment of theta(t) as a sinusoidal function at the top of the sheet. I included this to determine what a reasonable initial condition would be for thetadot(t), and I have used thetadot(0)=0.013 in the solve block. If the theta(t) defined at the top of the sheet is used throughout then this is not correct as the system would be forced continuously with this sinusoidal function. We really only need an inital velocity perturbation to get the system going to then see if unstable behaviour results.

The system response is expected to be a sinusoidal function that will vibrate at the natural frequency wo. The problem requires an initial disturbance thetadot(0) = some positive value, and the system will be stable/unstable depending on water height, damping etc.

If we elimiate this theta definition at the top of the sheet then F(t,s), pp(t,y), Fv(t), M(t) would all need to include theta as a variable.

When I included theta in F(t,s,theta) the symbolic solution gives an error as being recursive - which I guess makes sense. So I am bit stumped where to go from here?

Cheers

Ross

[Werner_E]

Yes, seems to be not that easy. Alan suggested to solve for F(), too and use the derivative of F in the solve block. But I am not quite sure, how.

[AlanStevens]

Werner Exinger wrote:

Yes, seems to be not that easy. Alan suggested to solve for F(), too and use the derivative of F in the solve block. But I am not quite sure, how.

If F(t,s) were not a function of s, then it would be straightforward (if a little tedious to develop!). Split F into the sum of F1 and F2 by expanding the cos(phio(t-tau)) term; differentate each of these terms wrt time to create functions for dF1/dt and dF2/dt which could be solved simultaneously with the dtheta/dt etc.

Unfortunately, F(t,s) is a function of s, so this doesn't quite work. I thought one possibility might be to use an average of F over s, but then noticed something strange - in the definition of pp the integrals over s go from zero to ul. However, for values of s less than mu^2/4g phio is imaginary! This doesn't feel right to me.

I guess I would need to look at the derivation of these equations in detail to have a chance of resolving this.

Alan

[remslie]

Alan,

If it helps, the paper where I have extracted the formulae is attached if you have the time/interest to examine. Equation 21 is the formula for the Dynamic Pressure pp(t,y). My formulation is slightly different as the paper assumes a flat gate shape and adjusts for curvature using an area coefficient cg in deriving the moment M(t). I have integrated the pressure over the curved surface to obtain the vertical force Fv(t) and then the moment M(t). This slight difference does not affect the the math approach. The integrals in s use an upper limit of infinity, whereas I have used ul = 100 - sensitivty seemed to indicate this value of ul was adequate.

I would be interested in your opinion if I have handled the evaluation of the last 2 integrals in equation 21 correctly. The last integral over dtau has the phase function cos(phio(s)), which is problematic as s is undefined in this integral. I evaluated this integral symbolically and then inserted the result in the second integral over ds.Is this "kosher" ?.

I agree with you that phio(s) is imaginary for small values of s. Being a phase function this is probably ok as it results in a change in the spatial pressure from positive to negative.

Regards

Ross

[AlanStevens]

I'll have a look. Might take me a couple of days though (and I can't guarantee I'll understand it!).

Alan

[remslie]

Alan,

Much obliged for your interest. I will soon be getting a chance to do some tests on the dam gate that exhibited unstable behaviour earlier this year hence my interest to try to simulate the problem.

Cheers

Ross

[remslie]

Alan,

I managed to find another paper that took a different approach to solving Equation 21. I will send my Mathcad sheet to you shortly out of interest to close out the problem.

I have a new problem with a new bug to fix! Vibration tests we have performed on a cardan shaft indicate it is exhibiting symptoms of a parametric instability. I am looking to solve several differential equations using a solve block. I have copied a previous odesolve method you provided some time ago but it is giving an errror I cannot seem to clear. I'm sure it will be trivial. If you can assist I would appreciate it.

Cheers

ross

[AlanStevens]

Works ok now - see attached.

Alan

PS Will be interested to see your other Mathcad sheet.

[remslie]

Alan,

Thanks again for your de-bugging expertise. I can't believe I missed the error with "t" in lieu of "tau". I mucked around for hours trying to source the problem!That's what a fresh pair of eyes does. With your modification I have the sheet working (attached) the method to plot results is clunky but it does the job.I will send through a copy of the Laplace transform sheet under a reply in the other discussion in case others haave an interest in the result.

Thanks again and have a great Xmas if we do not chat in the interim.

Ross

[AlanStevens]

Ross Emslie wrote:

... the method to plot results is clunky ...

There is no need to define K and G etc. Just transpose A and B in the return vector, for example. Also, no need for F at all; your definition of t works perfectly well - see attached.

Merry Christmas.

Alan

[remslie]

Alan,

There always seems to be one more bug to fix, although most of the time its me!

I have adapted your sheet to now include coupled 2nd order equations and it doesn't seeem to like it. The equations have been input as an expansion of the following matricies.It seems to solve OK when one of he second order terms is removed. from other discussions I have seen it appraes Mathcad can solve such coiupled systems, however, have I done something that is not allowed?

My sheet is attached.

Thanks ross

[AlanStevens]

Always best to get a single equation for each of the derivatives. You should check that I've rearranged things correctly. See attached for working version.

Alan

[remslie]

Alan,

That makes perfect sense when you see the manipulation you have done. Yes it checks out.

I can forge ahead now.

Thanks again for all your assistance this year. Have a great Xmas.

Cheers

ross