Yhtälöiden ratkaisemista Lambertin funktion avulla

Heikki Apiola
Aalto yliopisto, matematiikan ja systeemianalyysin laitos, lehtori emeritus
heikki.apiola@aalto.fi

Lähtökohta

Abitytöt Elle Oksanen ja Vilja Varis suorittivat vaihto-oppilasjaksoaan Utsjoen saamelaislukiossa maaliskuussa 2021. Sattuipa niin, että olivat ratkaisemassa GeoGebra-ohjelmalla yhtälöä, jossa esiintyi termejä \(x\) ja \(k^x\).

Kuvissa ovat tytöt pähkäilemässä tätä ja muitakin matematiikan ihmeitä vierailevan opettaja-Sissin avustuksella.

*ja vielä Wanhan Taskilan pihapiirissä.*

Kukaan ei muista, mikä tuo yhtälö tarkalleen ottaen oli, mutta katsotaanpa vaikka tätä: \(3^x = 2 x + 2\).

GeoGebra:n cas-laskin puhuu:

Mutta mikä ihmeen LambertW sinne ilmestyi, kysyivät tytöt.

Tässä kirjoituksessa valaistaan tuota salaperäistä LambertW-funktiota ja sen käyttömahdollisuuksia erityisesti yhtälöiden ratkaisemisessa.

Vähän historiaa

Johann Heinrich Lambertin elämä ajoittuu vuosille 1728–1777. Hänen kuuluisa aikalaisensa Leonhard Euler (1707–1783) työskenteli yhteistyössä Lambertin kanssa tämän tutkimuksen kohteena olevien yhtälöiden parissa etupäässä sarjamuotoisia ratkaisuja kehitellen.

Molemmat matemaatikot olivat ajan tyyliin erittäin laaja-alaisia ulottaen tutkimuksensa myös eri luonnontieteiden alueelle.

Lambertin tutkimusalueet liikkuivat mm. lukuteoriassa, tilastotieteessä, tähtitieteessä, optiikassa, filosofiassa. Hänet tunnetaan myös ensimmäisestä \(\pi\):n irrationaalisuuden todistuksesta, jota ovat sittemmin täsmentäneet ja yksinkertaistaneet monet matemaatikot aina viime aikoihin saakka. Lambertin matemaattisten julkaisujen joukosta löytyy myös tulevaisuuden visio koneesta, joka kykenisi suorittamaan matemaattista symbolien käsittelyä, siis verrattomasti kehittyneempää laitetta kuin Pascalin vuonna 1642 rakentama aritmetiikkakone.

Runsaat 200 vuotta Lambertin pähkäilyjä myöhemmin ilmestyivät ensimmäiset (‘muistiongelmaiset’) symbolialgebraohjelmistot, ja kas, nykyisin on ‘cas’ jokaisen koululaisenkin ulottuvilla.

Lambert tutki muotoa \(x = q + x^m\) olevaa yhtälöä kehittämällä ratkaisun sarjaksi. Euler kehitti sarjaratkaisun yhtälölle \(x e^x = a\), joka on muunnettavissa alkuperäisen Lambertin tutkiman yhtälön ratkaisuksi. Tämäpä on juuri se yhtälö, jonka ratkaisu määrittelee ‘Lambertin funktion’ siinä muodossa kuin sitä nykyisin käsitellään.

Kirjoituksessa käytettävät ohjelmistot

Tässä kirjoituksessa ei käsitellä Lambertin funktion arvojen laskemiseen tarvittavia numeerisia ratkaisumenetelmiä, vaan käytetään ohjelmistoja, joissa kyseinen funktio on valmiiksi ohjelmoituna.

Valtaosa laskuista suoritetaan Matlab/Octave:lla. Joitakin Maple-esimerkkejä on mukana symbolilaskennan alueella, ja pari näkymää myös GeoGebra-ratkaisuihin, kuten yllä.

Kirjoituksessa [HA1] on ohjelinkkien lisäksi perusteellinen alkuunpääsyohjeisto aina Octave-onlinen käyttöönottoa myöten.

Tämän kirjoituksen viitteissä on lisäksi hyvä Octave:n yleisopas [Oct].

Ohjelmakoodeja ja ajotuloksia, erityisesti niiden tuottamia kuvia on tekstissä mukana joiltakin osin hiukan riisuttuina. Täydelliset koodit kaikkine yksityiskohtineen on saatavissa viitteessä [HA2].

Käänteisfunktioita

Yhtälön ratkaiseminen yleisessä muodossa \(f(x) = y\) voidaan nähdä käänteisfunktion \(f^{-1}\) arvon laskemisena pisteessä \(y\). Tämä edellyttää käänteisfunktion tai sen haaran olemassaoloa pisteen \(y\) ympäristössä. Jos ratkaistavana on esimerkiksi yhtälö \(x^2 = 4,\) jolloin \(f(x) = x^2\), niin ratkaisut \(\pm 2\) sijaitsevat kaksihaaraisen käänteisfunktion \(\sqrt{y}\) positiivisella ja negatiivisella haaralla. Vaihdetaan \(y\):n ja \(x\):n roolit, jolloin kuvaaja muuttuu oikealle avautuvaksi paraabeliksi, ja positiivinen ja negatiivinen käänteisfunktion haara \(\pm \sqrt{x}\) näkyvät normaaliasentoisessa koordinaatistossa.

Logaritmi

Jos emme olisi kuulleet logaritmista, niin tämäkin GeoGebra-lasku olisi hämmentävä:

Johdatuksena pääaiheeseen onkin hyvä palauttaa mieleen logaritmin tarina, sitä tytötkin tarvitsivat, ennen kuin oppivat rakastamaan logaritmeja :-)

Funktio \(x\mapsto e^x\) on kasvava ja jatkuva, käänteisfunktio on siis olemassa. Minkähän nimen sille antaisimme? No, olisiko mitään, jos kutsuttaisiin sitä nimellä ‘luonnollinen logaritmi’? Kuulostaa monimutkaiselta, mutta olkoon menneeksi, otetaan sentään käyttöön helppo merkintä: ‘ln’. Siis jos \(y=e^x,\) niin \(x\) on \(y\):n luonnollinen logaritmi \(x=\ln(y)\).

Tässä siis pystyakselilla ovat lähtöarvot \(x\) ja vaaka-akselilla tulosarvot \(y=e^x\).

Kuvat syntyivät seuraavilla komennoilla Matlab:ssa/Octave:ssa:

x=linspace(-2,2,15);
y=exp(x);
plot(x,y); % Pisteet (x(k),y(k)).
hold on
plot(x,y,'*','MarkerSize',7);
title('y=e^x')
xlabel('x');ylabel('y')

Vaihdetaan x ja y plot-komennossa:

figure % Uusi grafiikkaruutu.
plot(y,x);  % Pisteet (y(k),x(k)).
grid on;hold on
plot(y,x,'*','MarkerSize',8);
title('x= ln y')
xlabel('y');ylabel('x')

Käänteisfunktion arvot saadaan siis valmiiksi lasketuissa funktion \(f(x) = e^x\) arvopisteissä \(y_k=f(x_k)\).

Käänteisfunktiokuvassa näkyvät merkityt (’*’) -pisteet ovat tarkkoja käänteisfunktion arvoja (\(f\)-arvojen laskennan tarkkuudella), koska ne ovat lukuja \(f(x_k)\) vastaavia, pystyakselille sijoitettuja lähtöarvoja \(x_k\).

Matlab/Octave suorittaa visuaalisesti paloittain lineaarisen interpolaation \((y_k,x_k)\)-pisteiden välillä. Jos halutaan laskea käänteisfunktion arvo mielivaltaisessa pisteessä \(y\), voidaan suorittaa numeerisesti paloittain lineaarinen (tai ‘palapolynomi’-) interpolaatio, jonka Matlab/Octave tekee visuaalisesti.

Menetelmää kutsutaan käänteiseksi interpolaatioksi ja sitä käytetään sopivin muunnelmin osana useita menetelmiä yhdistävissä yhtälöiden ratkaisemisen ‘hybridialgoritmeissa’. Tällä kertaa en etene pitemmälle numeeristen menetelmien parissa.

Kun funktio määritellään annetun funktion \(f\) käänteisfunktiona, on paikallaan miettiä, mitä johtopäätöksiä voidaan tehdä funktion \(f\) ominaisuuksien perusteella käänteisfunktion \(f^{-1}\) ominaisuuksista. Ainakin derivaatta voidaan johtaa.

Logaritmin derivaatta: Käänteisfunktion derivaatalle pätee \(D f^{-1}(y) = \frac{1}{f'(x)},\) missä \(x=f^{-1}(y)\). Jos \(f(x) = e^x\) niin \(f^{-1} = \ln\). Koska \(D\,e^x = e^x\), niin logaritmin derivaatta on \(D(\ln y) = \frac{1}{e^x},\) missä \(e^x = y\), joten \(D(\ln y) = \frac{1}{y}\).

Samalla mallilla voidaan johtaa kohta määriteltävän LambertW-funktion derivaatta. Muotoillaan harjoitustehtäväksi, kun sen aika koittaa.

Lambertin W-funktio

Toimitaan aivan samoin kuin edellä eksponenttifunktion ja logaritmin tapauksessa. Lähtökohtana on \(e^x\):n sijasta funktio

\[f(x)=x e^x.\]

Lähdetään siis etsimään käänteisfunktiota tälle. Piirretään aivan kuten edellä, nyt on \(e^x\):n sijasta \(x e^x\):

f= @(x) x.*exp(x) %Funktiomäärittely,huomaa(.*)
x=linspace(-6,1,1000);
y= f(x); 
plot(x,y,'k'); grid on
title('f(x)=xe^x')

Koska \(f'(x)= e^x(x+1),\) vahvistuu todeksi kuvan kertoma: \(f\):n minimi on kohdassa \(x=-1\), ja \(f\) on vähenevä vasemmalla ja kasvava oikealla puolella. Huomaa, että \(f(x)\to 0\), kun \(x\to -\infty\), ja \(f(x)\to \infty\), kun \(x\to \infty\). Minimin arvo \(f(-1) = -1/e \approx -0.367\). Siis jos \(-1/e < y < 0\), saadaan kaksi \(x\):n arvoa, joille \(x e^x = y\), ja jos \(y \geq 0\), saadaan yksi \(x\). Siis edellisillä \(y\):n arvoilla on kaksihaarainen käänteisfunktio ja jälkimmäisillä yksikäsitteinen.

Merkitään käänteisfunktiota tai sen haaroja \(W\):llä tai \(LambertW\):llä. Niinpä \(W\) suhtautuu funktioon \(x\mapsto xe^x\) aivan kuten luonnollinen logaritmi \(\ln\) suhtautuu eksponenttifunktioon \(x\mapsto e^x\).

Piirretään käänteisfunktio aivan kuten edellä ‘exp/log’-puuhassa. Vaihdetaan \(x\) ja \(y\), eli piirretään pisteet \((y_k,x_k)\). Jatketaan edellisen istunnon muuttujalla ‘x’.

x1=x(x>=-1);y1=f(x1);  %1. haara -> W_0
x2=x(x<=-1);y2=f(x2);  %2. haara -> W_{-1}
p1=plot(y1,x1,'b','LineWidth',2);
hold on;grid on
p2=plot(y2,x2,'r');
title('Käänteisfunktion W haarat')
legend([p1,p2],{'W_0','W_{-1}'})
e=exp(1); % e ei ole varattu symboli.
plot([-1/e -1/e],[-6,-1],'--k')
text(-1/e-0.1,-5.8,'-1/e')

Vaihdetaan käänteisfunktiossa y:n ja x:n roolit, ja kootaan yhteen edellä jo osittain todettua:

Käänteisfunktion kasvava haara \(W_0\) on määritelty joukossa \(x \ge -1/e\approx -0.3679\) ja saa arvoja \(y \ge -1\).
Pienenevä haara \(W_{-1}\) on määritelty välillä \(-1/e \leq x < 0\) ja se saa kaikki arvot \(y \leq -1\). Muista: \(\lim_{x\to -\infty} x e^x = 0\).
Välillä \(-1/e <x < 0\) käänteisfunktiolla on kaksi haaraa, vähenevä: \(y\leq -1\) ja kasvava: \(-1 \leq y < 0\).

Kompleksialueella \(W\) on äärettömän monihaarainen (kuten vaikkapa trigonometristen funktioiden käänteisfunktiot jo reaalialueella). Kompleksiset \(W\)-funktion haarat pohjautuvat kompleksialueen logaritmifunktion haaroihin, mutta niihin ei tämän kirjoituksen alue ulotu.

Merkinnällä \(W\) tarkoitetaan yleensä käänteisfunktion jotain haaraa, \(W_0\), \(W_{-1}\) tarkoittavat käänteisfunktion reaalista kasvavaa haaraa ja vähenevää haaraa vastaavasti. Useimmissa ohjelmistoissa \(W\)-funktiolla on nimenään ‘lambertW’, vaihtelevan kokoisin kirjaimin, ja ensimmäinen argumentti on reaalialueella haaraan viittaava \(0\) tai \(-1\). Jos 1. argumentti puuttuu, ts. annetaan vain yksi argumentti, niin käsitellään päähaaraa \(W_0\).

Nykyisin monet ohjelmistot ja kielet sisältävät Lambertin funktion koodin. Tässä on Octave-onlinen help-tekstin ote, jossa samalla kerrataan ohjelman merkinnöin edellä esitellyt funktion perusominaisuudet.

 octave:1> help lambertw
 -Function: lambertw(Z)
 -Function: lambertw(N, Z)
  This function satisfies W(z).*exp(W(z)) = z, 
  and can thus be used to express solutions of 
  equations involving exp's or log's.
  N must be integer, and specifies the branch of W.
  W(z) is a shorthand for W(0,z).  
  Branches 0 and -1 are the only ones that can 
  take on non-complex values.
  For example, the principal branch W(0,z) passes 
  through the point (0, 0):
    lambertw (0)   -> ans = 0
    lambertw(-1,0) -> ans = -Inf
  And the 0 and -1 branches coincide for the 
  real value:
     x = -1/exp(1);
     lambertw (0, x)   -> ans = -1
     lambertw (-1, x)  -> ans = -1

Esimerkkejä

Yhtälöt, jotka sisältävät termejä \(x^p\), \(e^x\), \(\ln x\), pyritään saattamaan muotoon \(v e^v = b\), jolloin \(W\)-funktion määritelmästä seuraa: \(v = W(b)\). (Vertaa tähän: Jos edellä kertoja \(v\) puuttuisi, niin ratkaisu olisi \(v = \ln b\).)

Esim. 1. Ratkaise alussa olleen GeoGebra-istunnon yhtälö \(3^x = 2 x +2\).

Ratkaisu: Yhtälössä esiintyy \(x\) ja \(3^x\), joten Lambertista voisi olla apua. Pyritään aluksi muotoon \(u 3^u\), missä \(u\) on sopiva apumuuttuja:

\[3^x = 2 x +2 \iff (x+1)3^{-x} = \frac{1}{2}.\]

Kerrotaan \(-\frac{1}{3}\):lla, jolloin \((-x-1)3^{-x-1} = -\frac{1}{6}.\)

Merkitään \(u=-x-1\), otetaan \(3\):n sijasta \(e\) kantaluvuksi ja kerrotaan yhtälö \(\ln 3\):lla, jotta päästään muotoon \((*)e^{(*)} = b\):

\[u(\ln 3) e^{u\ln 3} = - \frac{\ln 3}{6}.\]

Siis \(u \ln 3 = W(- \frac{\ln 3}{6})\).

Lasketaan \(-\frac{\ln 3}{6} = -0.1831\), joka kuuluu välille \((-\frac{1}{e}, 0)\), joten molemmat haarat \(W_0\), \(W_{-1}\) antavat reaalisen arvon.

Näin saadaan

\[\begin{aligned} x_{0} &=-1-\frac{1}{\ln\,3}W_0(\frac{-\ln 3}{6}),\\[4pt] x_{1} &=-1-\frac{1}{\ln\,3}W_{-1}(\frac{-\ln 3}{6}). \end{aligned}\]

Tämän tyyppisellä tekniikalla saadaan pienellä harjoittelurutiinilla yllättävän suuri joukko yhtälöitä, jotka sisältävät termejä \(x^p\), \(\ln x\), \(e^x\) ratkaistuksi Lambertin funktion avulla.

Lasketaan vielä likiarvot (Octave:ssa ‘log’ on luonnollinen logaritmi):

>> x0=-1-lambertw(0,-log(3)/6)/log(3)
   x0 = -0.7901  % Laskun tulos
>> x1=-1-lambertw(-1,-log(3)/6)/log(3)
   x1 = 1.4446   % Laskun tulos

Jatkotehtävä: cas-ratkaisu ja kuvasta tarkistaminen.

Katsotaan, miten Matlab:n symbolilaskenta selviää tehtävästä:

>> syms x
>> solve(3^x == 2*x + 2,x)
ans =
- lambertw(0, -log(3)/6)/log(3) - 1

Muuten oivallisesti, mutta toinen ratkaisu

- lambertw(-1, -log(3)/6)/log(3) - 1

jäi puuttumaan.

Osaisiko Maple?Annetaan LaTeX:n tulkata Maple:n antamat tulokset:

> ratk:=solve(3^x = 2*x + 2,x);
> latex(ratk[1]);

\[-\frac{W\! \left(-\frac{\ln \left(3\right)}{6}\right)+\ln \! \left(3\right)}{\ln \! \left(3\right)}\]

> latex(ratk[2]);

\[-\frac{W\! \left(-1, -\frac{\ln \! \left(3\right)}{6}\right)+\ln \! \left(3\right)}{\ln \! \left(3\right)}\]

Siispä Maple osasi ratkaista täydellisesti.

Kirjoituksen alkujohdannossa olevaa GeoGebra-kuvaa tihrustamalla näkyy, että GeoGebra:n cas pärjäsi tässä myös.

Esim. 2. Ratkaise yhtälö \(x=a e^x\). Millä \(a\):n arvoilla on 2, 1, 0 reaalista ratkaisua?

Ratkaisu: Kerrotaan yhtälö \(-e^{-x}\):llä, jolloin

\[-a = -x e^{-x}.\]

Lambertin \(W\)-funktion määritelmän mukaan saadaan \(x=-W(-a)\).

Katsotaan W-funktion kuvaa:

- 2 ratkaisua, kun \(-\frac{1}{e} \leq -a \leq 0\), eli \(0 \leq a \leq{\frac{1}{e}}\).
- 1 ratkaisu, kun \(-a > 0\), eli \(a < 0\).
- 0 reaalista ratkaisua, kun \(-a < -\frac{1}{e}\), eli \(a > \frac{1}{e}\).

Piirretään vielä kuva yhtälön molemmista puolista valitsemalla kultakin yllä mainitulta arvoalueelta näytteeksi \(a\):lle arvot 0.5, 0.2, \(-1\).

Kuva tehtiin näillä komennoilla:

>> x=linspace(-2,3);
>> plot(x,x,'--k')
>> hold on
% 1) a> 1/e = 0.3679 =>ei ratk.Esim: a=0.5
% 2) 0 < a < 1/e   => 2 ratk. Esim. a=0.2
% 3) a < 0   => 1 ratk. Esim. a=-1   
>> a=0.5; 
>> plot(x,a*exp(x),'b') % Ei ratk.
>> a=-1; 
>> plot(x,a*exp(x),'c','LineWidth',2)%1 rat
>> a=0.2;
>> plot(x,a*exp(x),'r','LineWidth',2)%2 rat
>> grid on;ylim([-5 5])  
>> legend('y=x','a=0.5 > 1/e',...
    'a=0.2','a=-1 < 0')
>> title('x=ae^x leikkauspisteet...')

Kuva tukee yllä olevia johtopäätöksiä.

Lasketaan seuraavaksi Octave:n ‘lambertw’-funktiolla nuo kuvassa näkyvät ratkaisut.

>> a=0.5
>>-lambertw(0,-a)  % ans = 0.7940 - 0.7701i
>>-lambertw(-1,-a) % ans = 0.7940 - 0.7701i

Ei reaalisia ratkaisuja.

>> a=0.2
>> x0=-lambertw(0,-a)  % ans = 0.2592
>> x1=-lambertw(-1,-a) % ans = 2.5426

2 reaalista ratkaisua.

>> a=-1
>> x0=-lambertw(0,-a) % ans = -0.5671
>> -lambertw(-1,-a) % ans = 1.5339 + 4.3752i

1 reaalinen ratkaisu.

Siispä kaikki on kuin kuvassa.

Analyyttinen vai numeerinen

Tässä on sopiva kohta palata samaan pohdiskeluun, jota harjoitettiin diffyhtälökirjoituksessa [HA1], kun ‘täydennysfunktioina’ olivat Airyn ja Besselin funktiot.

Tällä kertaa voitiin Lambertin \(W\)-funktion avulla määrittää parametri \(a\), jota tehtävässä kysyttiin. Pelkkien numeeristen yhtälöratkaisijoiden avulla oltaisiin jouduttu katselemaan numeerisen hakuammuskelun puita näkemättä analyyttisen ratkaisun tarjoamia metsiköitä. Toisaalta tehokkaat numeeriset menetelmät mahdollistavat Lambertin funktion hyötykäytön. Molempi siis parempi. Ja jälleen huomataan, että analyyttisen ratkaisun käsite laajenee luontevasti laajentamalla sopivasti ‘hyväksyttävien’ funktioiden joukkoa.

Esim. 3. Ratkaise yhtälö

\[x^y = y^x.\]

Tässä on kysymys implisiittisessä muodossa annetun yhtälön ratkaisemisesta muodossa \(y(x)\) tai \(x(y).\) Aika veikeä tehtävä, ehkä vaikeakin?

Ensinnäkin on selvää, että triviaaliratkaisu on \(y=x\), mutta mitä muuta? Kokonaislukuratkaisuja tuskin löytyy muita kuin \(2\) ja \(4\). Viimeinen toivo taitaa olla yritys turvautua ystäväämme Lambertiin. Miten tuossa lähdettäisiin liikkeelle? Seuraava osoittautuu kelvolliseksi.

Korotetaan yhtälö puolittain potenssiin \(\frac{1}{x}\):

\[x^{\frac{y}{x}}= y \implies e^{-\frac{\ln x}{x}y} = y^{-1}.\]

Kerrotaan puolittain \(y\):llä ja sitten järjestetään kertoimeksi sama kuin eksponentissa kertomalla \(-\frac{\ln x}{x}\):llä. Saadaan:

\[(-\frac{\ln x}{x}y) e^{-\frac{\ln x}{x}y} = -\frac{\ln x}{x}.\]

Näin saatiin vasen puoli muotoon \((*)e^{(*)}\), ja voidaan turvautua Lambertin apuun:

\[-\frac{\ln x}{x}y = W(-\frac{\ln x}{x}).\]

Siispä

\[y= - \frac{W(-\frac{\ln x}{x}) x}{\ln x}.\]

Miltä näyttää kuvaaja? Puolet kuvasta saadaan tähän tapaan:

x=linspace(.1,10);
y=-x./log(x).*lambertw(-log(x)./x);
plot(x,y)

Toinen puoli saadaan käyttämällä laskevaa haaraa \(W_{-1},\) jolloin täytyy rajoittua niihin \(x\):n arvoihin, joilla

\[-\frac{1}{e} < -\frac{\ln x}{x} < 0.\]

Vaihtoehtoisesti voitaisiin käyttää hyväksi symmetriaa, eli vaihtaa \(x\) ja \(y\), jolloin saataisiin kuvan puuttuva puolikas taas päähaaran \(W_0\) avulla. Kuva tehtiin edellisellä tyylillä (ja siihen jäi pikku virhe: \(W_1\) pitää lukea: \(W_{-1}\)).

Kuvasta (hiirikäden tarkkuudella) poimittujen arvojen avulla voit testata yhtälön toteutumista likimäärin.

Vertailun vuoksi voidaan piirtää myös implisiittimuodossa \(0\)-korkeuskäyrä ‘contour’-funktiolla, joka edustaa ‘raakaa voimaa’ verrattuna hienostuneeseen Lambert-käsittelyyn. Yksityiskohdat tästä, kuten kaikista muistakin kirjoituksen koodeista ovat viitteessä [HA2].

Mihin muuhun LambertW kelpaa?

Viite [CGHJK] on perusteellinen esitys aiheesta sisältäen suuren joukon matematiikan ja luonnontieteiden sovellusaloja. Mainitsen tässä vain yhden, joka liittyy sekä edelliseen esimerkkiin että Solmun numerossa \(1/2021\) olleeseen ‘Delfiinikirjoitukseen’ [HA1]. Siinä esittelin yhtenä ratkaisutekniikkana ‘muuttujien erottelun’. Mainitsin, että tämä menettely johtaa usein implisiittiseen (ratkaisemattomassa muodossa olevaan) yhtälöön \(y\):n ja \(x\):n välillä. Viitteessä [CGHJK] esitellään palamisilmiöön liittyvä ‘Model combustion problem’, jota kuvaa alkuarvotehtävä

\[y'=y^2(1-y), \ y(0)=\epsilon > 0.\]

Siinä näytetään muuttujien erottelulla aikaansaatu implisiittinen yhtälö, josta artikkelin kirjoittajat johtavat \(W\)-funktion avulla esitettävän eksplisiittisen ratkaisun.

Värikäs visuaalinen tiivistelmä Lambert-teemasta on nähtävissä ja seinätauluksi ladattavissa Lambert-posterissa:

http://www.orcca.on.ca/LambertW/

Käy ihmeessä katsomassa, se on ihan MUST!

Muutama sana W-funktion laskennasta

Julkaisu [CGHJK] on perusteellinen tietolähde koko aihepiiriin ‘Lambertin funktio’. Siinä kerrotaan, kuinka Maple-symbolilaskentaohjelmiston kehitysryhmä täydennettynä Donald Knuth’illa ryhtyi tutkimaan sen historiaa ja erilaisia käyttömahdollisuuksia sekä erityisesti kehittämään numeerista algoritmia, joka olisi kyllin tehokas ja tarkka, jotta se täyttäisi cas-ohjelman ‘mielivaltaisen tarkkuuden’ vaatimuksen ja toimisi tehokkaasti myös singulaaripisteen lähistöllä vieläpä funktion kompleksihaarojen suhteen. Julkaisussa pohditaan myös nimien ‘LambertW’ ja ‘W’ historiaa ja sen useitakin ‘uudelleen löytymisiä’ myös mm. \(\Omega\)-nimisenä. Matemaatikot Polya ja Szegö (1925) arvellaan myös nimen ‘W’ ensimmäisiksi käyttäjiksi.

Tekijät kehittivät ‘mielivaltaisen tarkkuuden’ laskenta-algoritmin, joka siis kirjoitettiin Maple-ohjelmistolle. Algoritmi on toteutettu ‘Halleyn menetelmällä’, joka on vähemmälle huomiolle jäänyt ‘terästys’ tunnetulle Newtonin menetelmälle. Julkaisu [CM] on helppolukuinen tiivistelmä, joka sisältää myös hyvin selkeän, Halleyn menetelmällä kirjoitetun Matlab-koodin ‘W’-funktiolle.

Nykyisin tämä funktio sisältyy moniin ohjelmointikieliin ja ohjelmistoihin nimellä ‘LambertW’ tai esim. ‘lambertw’, kuten ‘meidän’ Octave:ssa. cas-ohjelmistoista mainittakoon edellä mainittujen lisäksi Maxima, Mathematica ja Wolframin vastaava julkinen tuote Wolfram Alpha. (Kts. myös [WF].)

Derivaatta voidaan johtaa käänteisfunktion derivointikaavan avulla, kuten yllä logaritmin tapauksessa tai implisiittisen derivoinnin avulla kaavasta

\[W(x) e^{W(x)} = x.\]

Saadaan: \(W'(x)e^{W(x)}(1+W(x)) = 1\),

\[W'(x)=\frac{1}{e^{W(x)}(1+W(x))} = \frac{W(x)}{x(1+W(x))}.\]

Jälkimmäinen muoto saadaan sijoittamalla \(e^{W(x)}\):n paikalle \(\frac{x}{W(x)}\). Jos haluaisimme laskea arvon \(W'(0),\) pitäisi käyttää edellistä muotoa, josta saadaan \(W'(0)=1\). Kokeilemani cas-ohjelmat palauttavat jälkimmäisen muodon, jota ei voi käyttää arvolla \(x=0\). cas-ohjelmien puutteena on usein erikoispisteiden käsittely, mistä jokaisella ‘vakavalla’ käyttäjällä lienee kokemuksia. (Kts. myös [CJ].)

Tehtävä: Huomaa, että yllä \(W\) tarkoittaa kumpaa tahansa haaraa \(W_0\), \(W_{-1}.\) Piirrä derivaattafunktion kummankin haaran kuvaajat ja vertaa niitä edellä olevaan Lambertin funktion kuvaajaan. Valitse piirtoalue sopivasti väliä \(-\frac{1}{e} < x <0\) tapauksessa \(W_{-1}\) ja aluetta \(x > -\frac{1}{e}\) tapauksessa \(W_0\) reunoilta (reunalta) supistaen. Mallia voit ottaa tiedostoista [HA2].

Siis lopultakin: “Mikä ihmeen LambertW?”

Saivatkohan Elle ja Vilja valaistusta kysymykseen?

Viitteet

[CGHJK] Corless, Gonnet, Hare, Jeffrey, and D. E. Knuth: Advances in Computational Mathematics, Vol 5 (1996) 329–359, ladattavissa:
https://cs.uwaterloo.ca/research/tr/1993/03/W.pdf

[CM] Cleve Moler: https://blogs.mathworks.com/cleve/2013/09/02/the-lambert-w-function/

[Posteri] http://www.orcca.on.ca/LambertW/

[Wiki] https://en.wikipedia.org/wiki/Lambert_W_function

[WF] https://mathworld.wolfram.com/LambertW-Function.html

[HA1] https://matematiikkalehtisolmu.fi/2021/1/diffyhtaloita.pdf

[HA2] https://matematiikkalehtisolmu.fi/2021/2/skriptit/
Sisältää tiedostot ‘LambertSkriptit.m’ ja
‘LambertSkriptitLIVE.pdf’. Edellisessä on koodit (Matlab, Octave), jälkimmäinen Matlab:n tyylikkäällä ‘LiveEditorilla’ ajettu pdf.

[Oct] https://octave.org/doc/v6.2.0/

[CJ] R.M. Corless and D.J. Jeffrey,“Well, It Isn’t Quite That Simple”, SIGSAM Bulletin, 26, 3 (1992) 2–6.