python - algebraic constraint to terminate ODE integration with scipy -

-

i'm using scipy 14.0 solve system of ordinary differential equations describing dynamics of gas bubble rising vertically (in z direction) in standing still fluid because of buoyancy forces. in particular, have equation expressing rising velocity u function of bubble radius r, i.e. u=dz/dt=f(r), , 1 expressing radius variation function of r , u, i.e. dr/dt=f(r,u). rest appearing in code below material properties. i'd implement account physical constraint on z which, obviously, limited liquid height h. consequently implemented sort of z<=h constraint in order stop integration in advance if needed: used set_solout in order so. situation code runs , gives results, set_solout not working @ (it seems z_constraint never called actually...). know why? there more clever idea, may in order interrupt when z=h (i.e. final value problem) ? right way/tool or should reformulate problem?

thanks in advance

emi

from scipy.integrate import ode   db0 = 0.001 # init bubble radius y0, t0 = [ db0/2 , 0. ], 0. #init conditions h = 1  def y_(t,y,g,p0,rho_g,mi_g,sig_g,h):     r = y[0]      z = y[1]     z_ = ( r**2 * g * rho_g ) / ( 3*mi_g )  #velocity     r_ = ( r/3 * g * rho_g * z_ ) / ( p0 + rho_g*g*(h-z) + 4/3*sig_g/r ) #r dynamics         return [r_, z_]  def z_constraint(t,y):     h = 1  #should rather variable..     z = y[1]      if z >= h:         flag = -1      else:         flag = 0     return flag  r = ode( y_ ) r.set_integrator('dopri5') r.set_initial_value(y0, t0) r.set_f_params(g, 5*1e5, 2000, 40, 0.31, h) r.set_solout(z_constraint)  t1 = 6 dt = 0.1  while r.successful() , r.t < t1:     r.integrate(r.t+dt)

you're running this issue. set_solout work correctly, must called right after set_integrator, before set_initial_value. if introduce modification code (and set value g), integration terminate when z >= h, want.

to find exact time when bubble reached surface, can make change of variables after integration terminated solout , integrate respect z (rather t) z = h. paper describes technique m. henon, physica 5d, 412 (1982); may find this discussion helpful. here's simple example in time t such y(t) = 0.5 found, given dy/dt = -y:

import numpy np scipy.integrate import ode  def f(t, y):     """exponential decay: dy/dt = -y."""     return -y  def solout(t, y):     if y[0] < 0.5:         return -1     else:         return 0  y_initial = 1 t_initial = 0  r = ode(f).set_integrator('dopri5') r.set_solout(solout) r.set_initial_value(y_initial, t_initial)  # integrate until solout constraint violated r.integrate(2)  # new system t independent variable: see henon's paper details. def g(y, t):     return -1.0/y  r2 = ode(g).set_integrator('dopri5') r2.set_initial_value(r.t, r.y)  r2.integrate(0.5)  y_final = r2.t t_final = r2.y  # error: difference between found , analytical solution print t_final - np.log(2)

Comments

Post a Comment

Popular posts from this blog

C# random value from dictionary and tuple -

-
i have following code, , need random item contains string "region_1" dictionary. public static dictionary<int, tuple<string, int, currencytype>> itemarray = new dictionary<int, tuple<string, int, currencytype>>() { { 0xb3e, new tuple<string, int, currencytype>("region_1", 1500, currencytype.fame) } }; public static int generateitemid(string shopid) { var generateditem = itemarray.elementat(new random().next(0, itemarray.count)).key; } how select this? it's not possible efficiently... first create static random @ class level... prevent non-random behaviour if run query on short period of time... (it's seeded discrete clock) static random rnd = new random(); then: var item = itemarray.values .where(t => t.item1 == shopid) .orderby(_ => rnd.next()) .firstordefault()
Read more

cgi - How do I interpret URLs without extension as files rather than missing directories in nginx? -

-
so i'm trying install zoneminder , have run under nginx, has few compiled files need run using fast-cgi. these files lack extension. if file has extension, there no issue , nginx interprets file , return 404 if can't find it. if has no extension, assume it's directory , return 404 when can't find sort of index page. here have now: # security camera domain. server { listen 888; server_name mydomain.com; root /srv/http/zm; # enable php support. include php.conf; location / { index index.html index.htm index.php; } # enable cgi support. location ~ ^/cgi-bin/(.+)$ { alias /srv/cgi-bin/zm/; fastcgi_pass unix:/run/fcgiwrap.sock; fastcgi_param script_filename $document_root/$1; include fastcgi.conf; } } the idea is, under cgi-bin directory goes through fastcgi pipe. any idea on howi can fix this? upon closer inspection, realized zoneminder has 2 cgi sc...
Read more

.htaccess - htaccess convert request to clean url and add slash at the end of the url -

-
ok , searched more 4 days , more 20 posts here in stach or in websites without kind of success first here htaccess file rewriteengine on rewriterule ^([a-za-z0-9_]+)/?$ index.php?page=$1 [nc] rewriterule ^([a-za-z0-9_]+)\/([a-za-z0-9_]+)\/?$ index.php?page=$1&action=$2 [nc] i tried many codes 1 worked me can call request , works me www.example.com/pagename => www.example.com/index.php?page=pagename here questions , please provide me information "in detail" understand code well 1- how convert request automatically www.example.com/index.php?page=pagename www.example.com/pagename 2- if tried add slash @ end shows me page without images or styling mean www.example.com/pagename => works www.example.com/pagename/ => doesn't work that's :) for first part need 301 rule this: rewritecond %{the_request} /index\.php\?page=([^\s&]+) [nc] rewriterule ^ /%1? [r=301,l] for second part use absolute path in css, js, images files r...
Read more