I am trying to calculate the percentage of solar radiation that is in the visible spectrum.
My code is:
h=6.626e-34; %planck constant
c=3e8; %speed of light
k=1.38066e-23; %boltzman constant
T1=5777; %temperature in surface of sun
syms f
J1 =(2*pi*h/(c^2))*((f^3)/(exp((h*f)/(k*T1))-1)); %planck radiation law
hold on;
ezplot(J1,[0,1.5e15])
a=int(J1,f,0,inf) %total energy radiated
b=int(J1,f,4e14,8e14) %energy between 400-800Thz (visible radiation spectrum)
Vp = (b/a)*100 %percentage of visible/total radiation
As long as the graph is as it was expected, this means that I have not completely ruined it, the result of the integrals is as follows:
Vp = -1888568826285205004703258735621345426367059580820393216707788800000000000*((73396718487075910602267519716779133887030184268951416015625*log(1 - exp(23642358029674224853515625/7114894705749824515342336)))/205953985278202888163058116890711980917418253877248 - (73396718487075910602267519716779133887030184268951416015625*log(1 - exp(23642358029674224853515625/3557447352874912257671168)))/25744248159775361020382264611338997614677281734656 - (2390487322005187985890576650155251369405251302897930145263671875*polylog(2, exp(23642358029674224853515625/3557447352874912257671168)))/1857466834100924357302864708366291649120175241251782656 + (25952248522181378144831874533777514511468878135769288445192954296875*polylog(3, exp(23642358029674224853515625/3557447352874912257671168)))/67008813354583054015095330308295397033299967000474002915328 - (110058495767576691259417256590823904271799518678985866399923893187011219092083*polylog(4, exp(23642358029674224853515625/3557447352874912257671168)))/1888568826285205004703258735621345426367059580820393216707788800000000 + (2390487322005187985890576650155251369405251302897930145263671875*polylog(2, exp(23642358029674224853515625/7114894705749824515342336)))/7429867336403697429211458833465166596480700965007130624 - (25952248522181378144831874533777514511468878135769288445192954296875*polylog(3, exp(23642358029674224853515625/7114894705749824515342336)))/134017626709166108030190660616590794066599934000948005830656 + (110058495767576691259417256590823904271799518678985866399923893187011219092083*polylog(4, exp(23642358029674224853515625/7114894705749824515342336)))/1888568826285205004703258735621345426367059580820393216707788800000000 + 101409666798338314597227594049400067888200283050537109375/22835963083295358096932575511191922182123945984))/(12228721751952965695490806287869322696866613186553985155547099243001246565787*pi^4)
This is just a numeric expression (no constants), but I was expecting (and I'm trying to find) a single value.
Any ideas how to overcome this? thanks in advance
double(Vp)
returns 44.3072, which is exactly what I was looking for so far
vpa(Vp)
returns 44.307238264260285485868531074049 + 1.5008384323384489567473242679822e-35 * i, which is forwarded
norm(vpa(Vp))
giving the same result as double (Vp)
However, when I added a few more lines of code:
d=int(J1,f,0,4e14); %infrared energy
e=int(J1,f,8e14,inf); %ultraviolet energy
Ir = (d/a)*100; %percentage of infrared radiation
Uv = (e/a)*100; %percentage of ultraviolet radiation
double (Ir) :
Error using mupadmex
Error in MuPAD command: DOUBLE cannot convert the input expression into a double array.
If the input expression contains a symbolic variable, use the VPA function instead.
Error in sym/double (line 514)
Xstr = mupadmex('symobj::double', S.s, 0);
Error in symplanckradlaw (line 21)
Infrared=double(Ir)
vpa (Ir) (vpa (Ir)) ,
(abs(7.3340024900713941224341547032249e-56*limit((652255981594442743874468745505068648842285814001516463259648*f^2*polylog(2, exp((3873563939581825*f)/466281739436020499437475332096)))/15004497594028668398955870330625 - (608270107310811468411217054651916403272252179121058584901915330886243977872955782952124416*f*polylog(3, exp((3873563939581825*f)/466281739436020499437475332096)))/58120880811771703455030054556291506586990890625 - f^4/4 + (466281739436020499437475332096*f^3*log(1 - exp((3873563939581825*f)/466281739436020499437475332096)))/3873563939581825 + (283625243683820020974472587101002022201492492463545426687503953676410023626686775978867575742611811818507908158310055936*polylog(4, exp((3873563939581825*f)/466281739436020499437475332096)))/225134948049212098682315198853176286979563186266469146812890625, f == 0, Right) - 146.24549829953781858190522266202 - 2.501397387230748261245540446637e-36*i)^2)^(1/2)