Forums: » Register « |  Free Tools |  User CP |  Games |  Calendar |  Members |  FAQs |  Sitemap |  Support |

New Free Tools on Dev Shed!

#1
November 26th, 2012, 08:31 AM
 leonardo2887
Registered User

Join Date: Nov 2012
Posts: 2
Time spent in forums: 44 m 41 sec
Reputation Power: 0
[Solved] Problem with scipy.optimize.fsolve

Dear all,

I am new to this forum as well as new to python. I started writing a program which I wanted to use to fit some experimental data with a particular function. In order to write the function, a non linear equation has to be solved.
In the few code lines below I am trying just to print out the function together with the independent variable (Ctot). I will then change the print command with a curve_fit() commnad.

However, to the program:

python Code:
 Original - python Code
``` import scipyimport scipy.optimize, scipy.special, scipy.statsfrom math import sqrtfrom numpy import loadtxtfrom pylab import *from scipy.optimize import fsolve          def Sig(Ctot, K, u, DHNC, DHC):    CP = 0.5#Determines the free surfactant concentration from the binding isotherm and mass balance    def func(x, Ctot, CP, K, u):         S = x-0.5*(1+(K*u*(Ctot-CP*x)-1)/sqrt((1-K*u*(Ctot-CP*x))**2+4*K*(Ctot-CP*x)))         #print "S=", S         return S    Theta = fsolve(func, 0.5, args=(Ctot, CP, K, u))    print "Theta =", Theta    Cf = Ctot-CP*Theta    lamb = 1. + sqrt((K*Cf*Theta)/(1.-Theta)) #lambda0, eigenvalue from weight matrix    chi = (Cf*(1.-Theta)*K)/(Theta*lamb**2) #fraction of non-cooperatively bound surfactant    s = K*u*Cf #new variable        Theta1cf = 0.5*((u*K)/sqrt((s-1)**2+4*s/u)-((s-1)*(2*u*K*(s-1)+4*K))/(2*sqrt((s-1)**2+4*s/u)**3)) #first derivative of theta with respect to cf    a1 = -(4.0*Cf*Theta1cf*K)/(Theta*(sqrt((1-s)**2+4*s/u)+s+1)**2) #set of dummy functions to calculate S(Cf)    a2 = -(4.0*Cf*(1-Theta)*Theta1cf*K)/(Theta**2*(sqrt((1-s)**2+4*s/u)+s+1)**2)    a3 = +(4.0*(1-Theta)*K)/(Theta*(sqrt((1-s)**2+4*s/u)+s+1)**2)    a4 = -(8.0*Cf*(1-Theta)*K*((4*K-2*u*K*(1-s))/(2*sqrt((1-s)**2+4*s/u))+u*K))/(Theta*(sqrt((1-s)**2+4*s/u)+s+1)**3)    chi1cf = a1 + a2 + a3 + a4 #first derivative of chi with respect to Cf    DHm=DHNC*chi+(1-chi)*DHC #mean enthalpy of adsorption    return CP*(Theta1cf*DHm + Theta*(DHNC-DHC)*chi1cf)/(CP*Theta1cf+1)*0+Theta #calorimetric signal #   return Theta    K = 2.u = 10.DHNC = -1000DHC = -3500 Ctot, Y = loadtxt('JR400-SDS-05.dat', unpack=True)  print Ctot, Sig(Ctot, K, u, DHNC, DHC)
```

The equation which has to be solve is contained in func(). Exactly the same output is obntained if func()... is put outside the function Sig(). If I run the program I get the following output:
Code:
```[ 0.0203  0.028   0.0318  0.0355  0.0392  0.0429  0.0465  0.05    0.057
0.0604  0.0671  0.0775  0.1155  0.1532  0.1905  0.2273  0.2638  0.2999
0.3356  0.3709  0.4058  0.4404  0.4746  0.5084  0.5419  0.5751  0.6079
0.6404  0.6726  0.7044  0.7359  0.7671  0.798   0.8286  0.8589  0.8889
0.9186  0.948   0.9771  1.006   1.0346  1.0629  1.0909  1.1187  1.1462
1.2431  1.6506] Theta = [ 0.0233803]
[ 0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803]```

So, going with order, Ctot should be a ndarray. Also the output of Sig() is an ndarray with the same length as Ctot. What I really don't understand is, why is Theta an array of length 1?

I thought the problem may be, that Ctot is passed as an argument of the function, So I decided to print the values of the function (and uncommented the line print "S=", S in the function func() and the following output is optained:

Code:
```[ 0.0203  0.028   0.0318  0.0355  0.0392  0.0429  0.0465  0.05    0.057
0.0604  0.0671  0.0775  0.1155  0.1532  0.1905  0.2273  0.2638  0.2999
0.3356  0.3709  0.4058  0.4404  0.4746  0.5084  0.5419  0.5751  0.6079
0.6404  0.6726  0.7044  0.7359  0.7671  0.798   0.8286  0.8589  0.8889
0.9186  0.948   0.9771  1.006   1.0346  1.0629  1.0909  1.1187  1.1462
1.2431  1.6506] S= [ 0.5153606   0.51571418  0.51589425  0.51607325  0.51625595  0.51644247
0.51662769  0.51681139  0.51718988  0.5173792   0.51776316  0.51838895
0.52100977  0.52409536  0.52682602  0.52296529  0.45442974  0.00158271
-0.32610579 -0.41085581 -0.44222299 -0.45773801 -0.46679046 -0.4726678
-0.47678152 -0.47981323 -0.48212972 -0.48396008 -0.4854415  -0.48666077
-0.48768421 -0.48855512 -0.489305   -0.48995722 -0.49052955 -0.4910357
-0.49148641 -0.49189024 -0.49225406 -0.49258457 -0.4928851  -0.4931595
-0.49341099 -0.49364309 -0.49385718 -0.49450905 -0.49620455]
S= [ 0.5153606   0.51571418  0.51589425  0.51607325  0.51625595  0.51644247
0.51662769  0.51681139  0.51718988  0.5173792   0.51776316  0.51838895
0.52100977  0.52409536  0.52682602  0.52296529  0.45442974  0.00158271
-0.32610579 -0.41085581 -0.44222299 -0.45773801 -0.46679046 -0.4726678
-0.47678152 -0.47981323 -0.48212972 -0.48396008 -0.4854415  -0.48666077
-0.48768421 -0.48855512 -0.489305   -0.48995722 -0.49052955 -0.4910357
-0.49148641 -0.49189024 -0.49225406 -0.49258457 -0.4928851  -0.4931595
-0.49341099 -0.49364309 -0.49385718 -0.49450905 -0.49620455]
S= [ 0.5153606   0.51571418  0.51589425  0.51607325  0.51625595  0.51644247
0.51662769  0.51681139  0.51718988  0.5173792   0.51776316  0.51838895
0.52100977  0.52409536  0.52682602  0.52296529  0.45442974  0.00158271
-0.32610579 -0.41085581 -0.44222299 -0.45773801 -0.46679046 -0.4726678
-0.47678152 -0.47981323 -0.48212972 -0.48396008 -0.4854415  -0.48666077
-0.48768421 -0.48855512 -0.489305   -0.48995722 -0.49052955 -0.4910357
-0.49148641 -0.49189024 -0.49225406 -0.49258457 -0.4928851  -0.4931595
-0.49341099 -0.49364309 -0.49385718 -0.49450905 -0.49620455]
S= [ 0.51536061  0.51571419  0.51589426  0.51607325  0.51625596  0.51644248
0.5166277   0.5168114   0.51718988  0.51737921  0.51776317  0.51838895
0.52100978  0.52409537  0.52682603  0.52296529  0.45442977  0.00158278
-0.32610577 -0.4108558  -0.44222298 -0.457738   -0.46679045 -0.47266779
-0.47678151 -0.47981322 -0.48212972 -0.48396007 -0.48544149 -0.48666076
-0.4876842  -0.48855512 -0.48930499 -0.48995721 -0.49052954 -0.49103569
-0.4914864  -0.49189023 -0.49225405 -0.49258456 -0.49288509 -0.49315949
-0.49341099 -0.49364308 -0.49385717 -0.49450905 -0.49620454]
S= [-0.26435873 -0.38752518 -0.4518721  -0.51296691 -0.57023555 -0.62222561
-0.66713794 -0.70543147 -0.76762254 -0.79190657 -0.83083645 -0.87403956
-0.9478984  -0.97541116 -0.98909115 -0.9971163  -1.00237378 -1.00606038
-1.0087817  -1.01086954 -1.01252017 -1.01386032 -1.01496646 -1.01589448
-1.01668605 -1.017369   -1.0179624  -1.01848412 -1.01894631 -1.0193573
-1.01972617 -1.02005902 -1.02036081 -1.02063567 -1.020887   -1.02111766
-1.02133006 -1.02152626 -1.02170801 -1.02187739 -1.02203508 -1.02218222
-1.02231981 -1.02244918 -1.0225706  -1.02295327 -1.02405022]
S= [ 0.17850246  0.17856586  0.1784605   0.17824108  0.17787755  0.17733401
0.17658942  0.17561133  0.17264316  0.17054945  0.16462252  0.14785745
-0.17734815 -0.63439506 -0.75038486 -0.7880122  -0.80539998 -0.81515722
-0.82133763 -0.82558041 -0.82866368 -0.83100704 -0.8328412  -0.83431444
-0.8355262  -0.83653991 -0.83739768 -0.83813474 -0.8387747  -0.83933378
-0.83982772 -0.84026718 -0.84066063 -0.84101487 -0.84133542 -0.84162682
-0.84189282 -0.84213657 -0.8423607  -0.84256816 -0.84276006 -0.84293807
-0.84310361 -0.84325846 -0.8434031  -0.84385452 -0.8451126 ]
S= [ 0.10095452  0.09564024  0.09175446  0.08685868  0.08055503  0.07244711
0.06234627  0.04987657  0.01436871 -0.009486   -0.072742   -0.21701886
-0.6964262  -0.82038206 -0.85954789 -0.87728228 -0.88718851 -0.89343088
-0.89770125 -0.90079686 -0.90313922 -0.90497578 -0.9064492  -0.90765662
-0.90866632 -0.90952283 -0.91025624 -0.91089291 -0.91145066 -0.91194173
-0.91237859 -0.91276967 -0.91312175 -0.91344032 -0.9137299  -0.91399422
-0.91423643 -0.91445913 -0.91466456 -0.91485526 -0.91503215 -0.91519665
-0.91534998 -0.91549373 -0.91562827 -0.91604991 -0.917239  ]
S= [-0.16762416 -0.27002632 -0.33030406 -0.3923348  -0.45484414 -0.51514617
-0.56972298 -0.61780217 -0.69804371 -0.7297723  -0.78062969 -0.83639141
-0.92744069 -0.95915022 -0.97433616 -0.98304749 -0.98867162 -0.9925751
-0.99543487 -0.99761632 -0.99933314 -1.00072193 -1.00186475 -1.00282113
-1.00363513 -1.00433613 -1.00494423 -1.00547813 -1.00595052 -1.00637012
-1.00674634 -1.00708551 -1.0073928  -1.00767245 -1.00792799 -1.00816237
-1.00837808 -1.00857722 -1.00876161 -1.00893338 -1.00909322 -1.0092423
-1.00938167 -1.00951266 -1.00963557 -1.01002265 -1.01113017]
S= [  4.76482953e-02   2.97678644e-02   1.68644513e-02   7.33096042e-04
-1.98074383e-02  -4.57380447e-02  -7.70691621e-02  -1.13990120e-01
-2.07421266e-01  -2.60821382e-01  -3.73299783e-01  -5.35616150e-01
-8.10004282e-01  -8.78532934e-01  -9.04728036e-01  -9.17938564e-01
-9.25805375e-01  -9.30973544e-01  -9.34613166e-01  -9.37308281e-01
-9.39380972e-01  -9.41026981e-01  -9.42361212e-01  -9.43463877e-01
-9.44392526e-01  -9.45185056e-01  -9.45867199e-01  -9.46462027e-01
-9.46985175e-01  -9.47447387e-01  -9.47859846e-01  -9.48230106e-01
-9.48564267e-01  -9.48867310e-01  -9.49143340e-01  -9.49395772e-01
-9.49627471e-01  -9.49840852e-01  -9.50037973e-01  -9.50221208e-01
-9.50391382e-01  -9.50549816e-01  -9.50697661e-01  -9.50836403e-01
-9.50966388e-01  -9.51374516e-01  -9.52531912e-01]
S= [-0.03356636 -0.08413044 -0.11903666 -0.16038587 -0.20910712 -0.26443916
-0.3229645  -0.38206802 -0.49734822 -0.54827914 -0.63453103 -0.7317376
-0.87962475 -0.9236086  -0.94283128 -0.95328665 -0.95981079 -0.96423383
-0.96741938 -0.96981812 -0.97168695 -0.97318646 -0.97441219 -0.97543227
-0.97629641 -0.97703759 -0.97767831 -0.97823911 -0.97873398 -0.97917249
-0.97956482 -0.97991784 -0.98023712 -0.98052722 -0.98079192 -0.98103439
-0.98125727 -0.9814628  -0.98165291 -0.98182983 -0.98199432 -0.98214762
-0.9822908  -0.98242529 -0.98255139 -0.98294799 -0.98407812]
S= [ 0.00642834 -0.0268322  -0.05048608 -0.07947446 -0.1152534  -0.15839015
-0.2073521  -0.26070825 -0.37731261 -0.43430441 -0.53793338 -0.66269849
-0.85319019 -0.9055521  -0.92728945 -0.93877174 -0.94580673 -0.95051724
-0.95387972 -0.9563948  -0.95834411 -0.95990172 -0.96117065 -0.9622237
-0.96311367 -0.96387546 -0.96453282 -0.96510733 -0.9656136  -0.96606168
-0.96646215 -0.96682215 -0.96714745 -0.9674428  -0.9677121  -0.96795861
-0.96818507 -0.9683938  -0.96858676 -0.96876625 -0.96893306 -0.96908845
-0.96923353 -0.96936975 -0.96949744 -0.96989874 -0.97104002]
S= [  7.48916438e-04  -3.48837053e-02  -6.01406021e-02  -9.09675784e-02
-1.28792603e-01  -1.74032082e-01  -2.24868399e-01  -2.79634740e-01
-3.97108288e-01  -4.53518931e-01  -5.54699510e-01  -6.74836917e-01
-8.57629146e-01  -9.08503044e-01  -9.29805526e-01  -9.41113131e-01
-9.48062147e-01  -9.52724627e-01  -9.56057706e-01  -9.58553531e-01
-9.60489564e-01  -9.62037612e-01  -9.63299446e-01  -9.64347093e-01
-9.65232834e-01  -9.65991251e-01  -9.66645898e-01  -9.67218170e-01
-9.67722582e-01  -9.68169102e-01  -9.68568246e-01  -9.68927107e-01
-9.69251432e-01  -9.69545928e-01  -9.69814482e-01  -9.70060336e-01
-9.70286215e-01  -9.70494422e-01  -9.70686922e-01  -9.70865997e-01
-9.71032426e-01  -9.71187476e-01  -9.71332251e-01  -9.71468192e-01
-9.71595621e-01  -9.71996159e-01  -9.73135638e-01]
S= [ -1.90751672e-05  -3.59756293e-02  -6.14500778e-02  -9.25247598e-02
-1.30622514e-01  -1.76137918e-01  -2.27214523e-01  -2.82154814e-01
-3.99711351e-01  -4.56032066e-01  -5.56875738e-01  -6.76405949e-01
-8.58208639e-01  -9.08890700e-01  -9.30136771e-01  -9.41421630e-01
-9.48359424e-01  -9.53015625e-01  -9.56344855e-01  -9.58838157e-01
-9.60772449e-01  -9.62319243e-01  -9.63580146e-01  -9.64627083e-01
-9.65512268e-01  -9.66270243e-01  -9.66924532e-01  -9.67496511e-01
-9.68000678e-01  -9.68446993e-01  -9.68845962e-01  -9.69204674e-01
-9.69528871e-01  -9.69823254e-01  -9.70091710e-01  -9.70337478e-01
-9.70563280e-01  -9.70771420e-01  -9.70963858e-01  -9.71142878e-01
-9.71309258e-01  -9.71464263e-01  -9.71608998e-01  -9.71744902e-01
-9.71872298e-01  -9.72272735e-01  -9.73411976e-01]
S= [  5.51917419e-08  -3.59484214e-02  -6.14174490e-02  -9.24859642e-02
-1.30576937e-01  -1.76085493e-01  -2.27156151e-01  -2.82092157e-01
-3.99646722e-01  -4.55969707e-01  -5.56821784e-01  -6.76367066e-01
-8.58194262e-01  -9.08881076e-01  -9.30128545e-01  -9.41413968e-01
-9.48352040e-01  -9.53008397e-01  -9.56337723e-01  -9.58831087e-01
-9.60765423e-01  -9.62312248e-01  -9.63573174e-01  -9.64620129e-01
-9.65505328e-01  -9.66263313e-01  -9.66917612e-01  -9.67489598e-01
-9.67993770e-01  -9.68440091e-01  -9.68839065e-01  -9.69197780e-01
-9.69521980e-01  -9.69816366e-01  -9.70084824e-01  -9.70330595e-01
-9.70556399e-01  -9.70764540e-01  -9.70956980e-01  -9.71136001e-01
-9.71302382e-01  -9.71457389e-01  -9.71602124e-01  -9.71738029e-01
-9.71865426e-01  -9.72265865e-01  -9.73405113e-01]
S= [  4.05550246e-12  -3.59484999e-02  -6.14175431e-02  -9.24860761e-02
-1.30577069e-01  -1.76085644e-01  -2.27156319e-01  -2.82092338e-01
-3.99646908e-01  -4.55969887e-01  -5.56821940e-01  -6.76367178e-01
-8.58194303e-01  -9.08881104e-01  -9.30128568e-01  -9.41413990e-01
-9.48352062e-01  -9.53008418e-01  -9.56337743e-01  -9.58831108e-01
-9.60765443e-01  -9.62312268e-01  -9.63573194e-01  -9.64620149e-01
-9.65505348e-01  -9.66263333e-01  -9.66917632e-01  -9.67489618e-01
-9.67993790e-01  -9.68440111e-01  -9.68839085e-01  -9.69197800e-01
-9.69522000e-01  -9.69816386e-01  -9.70084844e-01  -9.70330615e-01
-9.70556418e-01  -9.70764559e-01  -9.70957000e-01  -9.71136021e-01
-9.71302402e-01  -9.71457409e-01  -9.71602144e-01  -9.71738049e-01
-9.71865445e-01  -9.72265885e-01  -9.73405132e-01]
S= [ -3.46944695e-18  -3.59484999e-02  -6.14175431e-02  -9.24860761e-02
-1.30577069e-01  -1.76085644e-01  -2.27156319e-01  -2.82092338e-01
-3.99646908e-01  -4.55969887e-01  -5.56821940e-01  -6.76367178e-01
-8.58194303e-01  -9.08881104e-01  -9.30128568e-01  -9.41413990e-01
-9.48352062e-01  -9.53008418e-01  -9.56337743e-01  -9.58831108e-01
-9.60765443e-01  -9.62312268e-01  -9.63573194e-01  -9.64620149e-01
-9.65505348e-01  -9.66263333e-01  -9.66917632e-01  -9.67489618e-01
-9.67993790e-01  -9.68440111e-01  -9.68839085e-01  -9.69197800e-01
-9.69522000e-01  -9.69816386e-01  -9.70084844e-01  -9.70330615e-01
-9.70556418e-01  -9.70764559e-01  -9.70957000e-01  -9.71136021e-01
-9.71302402e-01  -9.71457409e-01  -9.71602144e-01  -9.71738049e-01
-9.71865445e-01  -9.72265885e-01  -9.73405132e-01]
Theta = [ 0.0233803]
[ 0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803  0.0233803
0.0233803  0.0233803  0.0233803  0.0233803  0.0233803]```

So, also S, which is calculated in func() is an array of the same length of Ctot. What do I do wrong? Why can't I get Theta as an array with the same length of Ctot, with the value of the n-th element resulting from the n-th value of Ctot?

I am sorry if these questions are trivial, but I can't really understand what I am doing wrong.

With best regards,
Leonardo

#2
November 26th, 2012, 09:57 AM
 leonardo2887
Registered User

Join Date: Nov 2012
Posts: 2
Time spent in forums: 44 m 41 sec
Reputation Power: 0
Problem solved with by substituting line 19

Theta = fsolve(func, 0.5, args=(Ctot, CP, K, u))

with

Theta = fsolve(func, Ctot/2., args=(Ctot, CP, K, u))

The size of theta is given by the size of the initial guess, and as I put only one number, the size of theta was only 1.

With best regards,
Leonardo
b49P23TIvg agrees: Self solve. Yeah!!!

 Viewing: Dev Shed Forums > Programming Languages > Python Programming > Problem with scipy.optimize.fsolve