Page 1 sur 1
fsolve
Publié : sam. févr. 27, 2010 7:45 am
par Guillaume
Bonjour,
j'ai un petit problème avec fsolve que je ne rencontrais pas avant :
g(x):=-exp(2-x)-x*(-(exp(2-x)))+1/x
fsolve(g(x),x,2.1,newton_solver)
renvoie une approximation de l'une ou l'autre des 2 solutions ou undef ou 5.90931518923+4.54385779848*i
ou
terminate called after throwing an instance of 'std::bad_alloc'
what(): std::bad_alloc
Abandon
Re: fsolve
Publié : sam. févr. 27, 2010 7:53 am
par parisse
je n'ai pas de plantage, j'obtiens des racines complexes, ou 0.58... ou undef. C'est normal que newton ne renvoie pas toujours la même chose, étant donné le point de départ, la 1ère recherche échoue et la suivante utilise le générateur aléatoire pour déterminer un nouveau point de départ. Pour pouvoir reproduire l'abandon, il faudrait connaitre la valeur de seed du generateur aleatoire (il suffit de faire rand() juste avant la commande fsolve pour le connaitre).
Re: fsolve
Publié : sam. févr. 27, 2010 8:02 am
par Guillaume
Il a fallu une quinzaine d'essais :
Code : Tout sélectionner
31>> rand()
2108070049
// Time 0
32>> fsolve(g(x),x,2.1,newton_solver)
terminate called after throwing an instance of 'std::bad_alloc'
what(): std::bad_alloc
Abandon
Re: fsolve
Publié : sam. févr. 27, 2010 2:34 pm
par parisse
je viens d'essayer
for j from 1 to 15 do rand() od
rand() me renvoie alors 2108070049
puis g(x):=-exp(2-x)-x*(-(exp(2-x)))+1/x
fsolve(g(x),x,2.1,newton_solver) me renvoie alors 5.90931518923+4.54385779848*i
Vous avez une config spéciale par rapport à la config par défaut?
Re: fsolve
Publié : sam. févr. 27, 2010 3:02 pm
par Guillaume
C'est bizarre de tomber exactement sur la même "graine", non ?
Je n'ai pas de config spéciale. Un autre utilisateur m'avait signalé ce problème.
Re: fsolve
Publié : sam. févr. 27, 2010 3:33 pm
par parisse
non, c'est normal car je ne fais pas d'initialisation du générateur aléatoire à l'ouverture de session (on peut bien sur le faire explicitement avec srand(NULL) par exemple).
Re: fsolve
Publié : dim. févr. 28, 2010 11:33 am
par Guillaume
Bonjour,
je remplace newton par steffenson mais de nouveaux problèmes arrivent avec for ... in ... do ...end_for
Code : Tout sélectionner
TVapp(L,F,f):={
local k,j,S,m,fp;
f:=unapply(f,x);
fp:=function_diff(f);
S:=NULL;
j:=[seq(k,k=floor(L[0])..100)] minus F;
for k in j do for(m:=-5;m<=5;m++){ S:=S,resoudre_numerique(fp(x),x,k+m*0.1,steffenson_solver)};end_for
return(S)
}:;
Alors TVapp([0,+infinity],[0],ln(x))
ne fonctionne pas.
Quand je fais ligne par ligne avec giac ça marche
Code : Tout sélectionner
26>> L:=[0,+infinity];F:=[0];S:=NULL;
[0,+(infinity)],[0]
// Time 0
27>> j:=[seq(k,k=floor(L[0])..100)] minus F:;
"Done"
// Time 0
28>> for k in j do for(m:=-5;m<=5;m++){ S:=S,resoudre_numerique(ln(x),x,1+m*0.1,steffenson_solver)};end_for;
Re: fsolve
Publié : dim. févr. 28, 2010 2:51 pm
par parisse
je suppose que l'appel à la GSL est abusé par le changement de signe en x=0. Je pense qu'il vaudrait mieux utiliser une méthode à priori moins efficace de résolution, par exemple dichotomie (bisection_solver) si un changement de signe est repéré entre 2 points de la discrétisation. Les méthodes efficaces ont en effet un inconvénient c'est qu'il faut déjà être assez prêt d'une solution pour que ça converge (et ça converge alors très vite).
Re: fsolve
Publié : dim. févr. 28, 2010 6:29 pm
par Guillaume
Ce qui me trouble c'est que la même action mise dans une procédure bloque alors qu'elle fonctionne en dehors (dans l'une, un bad argument value et pas dans l'autre) : cela semble indépendant de la méthode de résolution.
Re: fsolve
Publié : dim. févr. 28, 2010 7:42 pm
par parisse
oui, bizarre... Je n'ai pas observé le bad arg value (mais je n'ai pas testé avec la 0.8.5, seulement avec la 0.9.0)
Re: fsolve
Publié : dim. févr. 28, 2010 8:14 pm
par Guillaume
N'y a-t-il pas de version 0.9 disponible ?....
Re: fsolve
Publié : lun. mars 01, 2010 8:05 am
par parisse
Seulement en source. Comme il contient une implementation de l'algorithme de PGCD modulaire pour les extensions algebriques, j'attend que le code soit bien stabilise avant d'en faire des versions binaires.
Re: fsolve
Publié : lun. mars 01, 2010 8:07 am
par parisse
je viens de tester avec la 0.8.5, et j'obtiens en reponse une liste de 0.0. Il y a donc bien quelque chose de mysterieux!