//////////////////////////////////////////////////////// // Méthode des lignes de niveaux pour la minimisation // // de la compliance d'une console // // // // Ecole Polytechnique, MAP 562 // // Copyright G. Allaire, C. Dousset 2007 // //////////////////////////////////////////////////////// // La fonction ligne de niveaux est appelée phi. La // forme est caractérisée par phi<0. La fonction phi // est convectée par un schéma caractéristique implicite // (cf. J. Strain, JCP 151, 1999). //////////////////////////////////////////////////////// int n=2; //paramètre de taille du maillage real mu=8; //coefficient de Lamé real lambda=1; //coefficient de Lamé real g1=0; //composante horizontale de la force appliquée sur le bord real g2=-5; //composante verticale de la force appliquée sur le bord real lag=0.5; //multiplicateur de Lagrange pour le poids real T=.1; //pas de descente ou pas de temps de la convection real TT=0.01; //pas de temps pour la réinitialisation real weak=0.0001; //coefficient multiplicatif du matériau "mou" représentant le vide real eps1=0.001; //petit paramètre pour ne pas diviser par zéro dans le calcul de la normale real eps2=0.001; //paramètre de régularisation de la vitesse real tol=0.0001; //tolérance pour la remontée de la fonction objectif real compliance; real poids; real objectif; real objectifinitial; real compliancem; real poidsm; real objectifm; int kx=4; //nombre de trous dans la direction horizontale de la forme initiale int ky=3; //nombre de trous dans la direction verticale de la forme initiale int N=30; //nombre d'itérations de la méthode de gradient int Ninit=5; //nombre d'itérations de la réinitialisation int iter; int iterinit; string legende; string sauve="console-ls"; ofstream object("console-ls.obj"); //Définition des bords et du maillage // bord droit du maillage border bd1(t=-0.5,-0.05) { x=1; y=t; label=1; }; border bd2(t=-0.05,0.05) { x=1; y=t; label=3; }; border bd3(t=0.05,0.5) { x=1; y=t; label=1; }; // bord gauche du maillage border bg(t=0.5,-0.5) { x=-1; y=t; label=2; }; // bord supérieur du maillage border bs(t=1,-1) { x=t; y=0.5; label=1; }; // bord inférieur du maillage border bi(t=-1,1) { x=t; y=-0.5; label=1; }; // signification des labels: 1=Neumann, 3=force appliquée, 2=Dirichlet mesh Th= buildmesh (bd1(9*n)+bd2(2*n)+bd3(9*n)+ bg(20*n)+bs(40*n)+bi(40*n)); //Définition des espaces d'éléments finis fespace Vh(Th,[P1,P1]); fespace Vh1(Th,P1); Vh [u,v],[w,s],[um,vm]; Vh1 N1, N2, phi, V, phim, M1, M2, S, phiinit, nabla, Xapprox, X, Vreg, vv, Xapproxm, Xm; // Définition de la forme initiale (perforée par des trous selon kx et ky) func phi0=-0.1+sin(pi*kx*x)*sin(pi*ky*(y-0.5)); phi=phi0; plot(phi,wait=0,fill=1,value=1,ps=sauve+".init.eps"); // Réinitialisation de phi pour que phi soit proche de la fonction distance // signée au bord. Résolution implicite linéarisée de : // (d phi)/(d t) + sign(phi) ( |gradphi)| - 1 ) = 0 Vh1 h1=hTriangle; real h=h1[].max; //taille maximale des triangles du maillage cout<<"taille maximale des mailles h="<objectif*(1+tol)) { T=T/2; cout <<"refus du pas de descente"<