!************************************************************************** !* Solve Laplace Equation by relaxation Method: d2T/dx2 + d2T/dy2 = 0 (3) * !* ---------------------------------------------------------------------- * !* Example #3: Determine the temperatures T(x,y) in a domain [x,y] made * !* of one square of length C with the following conditions: * !* along OC: T=0, along AB: T=1, along OA and CB: dT/dy = 0 * !* (see sketch below). * !* Main parameters: * !* C: length of square domain. * !* Number of subdivisions 8 in ox and 8 in Oy. * !* integration steps dx=C/8. * !* weight coefficient w = 1.3 * !* Maximum residual error ER = 0.001 * !* Maximum number of iterations IT = 100 * !* ---------------------------------------------------------------------- * !* SAMPLE RUN: * !* * !* Niter Max. Residual W * !* 43 0.000860 1.30 * !* Temperature * !* 0.000 0.125 0.249 0.374 0.499 0.624 0.749 0.875 1.000 * !* 0.000 0.125 0.249 0.374 0.499 0.624 0.749 0.875 1.000 * !* 0.000 0.124 0.249 0.374 0.499 0.624 0.749 0.875 1.000 * !* 0.000 0.124 0.249 0.374 0.499 0.624 0.749 0.875 1.000 * !* 0.000 0.124 0.249 0.374 0.499 0.624 0.749 0.875 1.000 * !* 0.000 0.124 0.249 0.374 0.499 0.624 0.749 0.875 1.000 * !* 0.000 0.124 0.249 0.374 0.499 0.624 0.749 0.875 1.000 * !* 0.000 0.124 0.249 0.373 0.498 0.624 0.749 0.875 1.000 * !* 0.000 0.124 0.249 0.373 0.498 0.624 0.749 0.874 1.000 * !* * !* ---------------------------------------------------------------------- * !* REFERENCE: "Méthode de calcul numérique- Tome 2 - Programmes en Basic * !* et en Pascal By Claude Nowakowski, Edition du P.S.I., * !* 1984" [BIBLI 04]. * !* * !* F90 Release By J-P Moreau, Paris. * !* (www.jpmoreau.fr) * !************************************************************************** ! Note: Here, what is new is the limit conditions involving partial derivatives dT/dy. ! Inside domain, the classical scheme is applied to calculate the residual value at current ! point (i,j): 4*T(i,j)-T(i-1,j)-T(i+1,j)-Ti,j-1)-T(i,j+1). Must be near zero (< ER). ! For limit lines OA and CB, the partial derivative is approximated by: ! dT/dy = (T(i,j+1)-T(i,j-1))/(2*dy) ! But the point T(i,j+1) for segment CB is out of the domain and becomes a fictitious point, ! so we apply the following scheme: r=4*T(i,j)-2*T(i,j-1)-T(i-1,j)-T(i+1,j). ! For the other segment OA, it becomes: r=4*T(i,j)-2*T(i,j+1)-T(i-1,j)-T(i+1,j). ! The exact solution for this problem is T=K*x (0