real pi=4*atan(1.); real alpha = 7*pi/4; real omega = pi/alpha; int N=8; // Le premier pas du maillage est h=1/N int Nit = 5; // Nbre de divisions par 2 du pas de maillage (a augmenter au maximum des possibilites de l'ordinateur) //Definition de la frontiere du domaine, label=0 correspond a la condition de Neumann, label=1 a la condition de Dirichlet. border D1(t=0,1){x=t; y=0; label=1;}; border D2(t=1,0){x=t*cos(alpha); y=t*sin(alpha); label=1;}; border C(t=0 , alpha){x=cos(t); y=sin(t); label=0;}; real[int] errH10(Nit), errL2(Nit), ordredeconvergence(Nit - 1), tabh(Nit); // tableaux d'erreur ,d'ordre de convergence, de pas de maillage // Definitions de fonctions utiles : r,theta (coordonnees polaires), phi et ses derivees, u0 et ses derivees et f. func r=sqrt(x^2 + y^2) ; func phi = (r<0.5) + (r>=0.5)*0.5*(1 + cos(2*pi*(r - 0.5))); func phiprime = -pi*(r>=0.5)*sin(2*pi*(r - 0.5)); func phiseconde = (-2*pi^2)*(r>=0.5)*cos(2*pi*(r - 0.5)); func rnot0 = r+exp(-100000*r); // Pour eviter les divisions par 0 en r=0 dans le calcul de theta. func theta = acos(x/rnot0) * (y>=0) + (2*pi - acos(x/rnot0))*(y<0); func u0 = phi*r^omega * sin(omega*theta); func GRADru0 = r^(omega - 1) *(phiprime*r + omega*phi)*sin(omega*theta); func GRADthetau0 = (omega/rnot0)*phi*r^omega*cos(omega*theta); func dxu0 = cos(theta)*GRADru0 - sin(theta)*GRADthetau0; func dyu0 = sin(theta)*GRADru0 +cos(theta)*GRADthetau0; func f = -r^(omega - 1)*((1 + 2*omega)*phiprime + r*phiseconde)*sin(omega*theta); ////////////////////////////// BOUCLE SUR LE PAS DE MAILLAGE /////// for(int i=0; i