Browse code

Исправление работы функций GSLL: преобразование в векторы возвращаемых значений

Dmitrii Kashin authored on 11/05/2015 15:15:10
Showing 6 changed files
... ...
@@ -26,7 +26,7 @@
26 26
   ((dimensions :initarg :dimensions
27 27
 	       :initform 2
28 28
 	       :accessor dimensions)
29
-   (headers/cell2d-xv :initform '()
29
+   (headers/cell2d-cv :initform '()
30 30
 		      :allocation :class)))
31 31
 
32 32
 (print-functions-for-class cell2d-cv 
... ...
@@ -78,6 +78,18 @@
78 78
 
79 79
 ;; вспомогательные функции
80 80
 
81
+(defun vector-double-float->vector (agrid)
82
+  "В документации сказано (см. раздел 11 системы antik), что объекты
83
+типа `vector-double-float' можно преобразовать к типу `vector'
84
+автоматически при помощи функции `cl-array' только в том случае, если
85
+эти объекты не были порождены cffi-функцией. К сожалению, у меня
86
+именно такой случай, поэтому есть необходимость разбирать поля
87
+полученных структур вручную."
88
+  (list->vector
89
+   (loop for i to (1- *number-of-conservative-variables*)
90
+	       for v = (cons (grid:aref agrid i) v)
91
+	       finally (return (reverse v)))))
92
+
81 93
 (defun curryr (function &rest args)
82 94
   (lambda (x) (apply function x args)))
83 95
 
... ...
@@ -20,7 +20,7 @@
20 20
 
21 21
   (let ((stop-flag nil)
22 22
 	(newgrid (clone grid)))
23
-
23
+    
24 24
     ;; сначала экспортируем начальное поле
25 25
     (export-grid grid :print-out-border print-out-border)
26 26
     ;; если finish-step установлен в 0, то больше ничего не делаем
... ...
@@ -60,13 +60,17 @@
60 60
 	    (when (>= (1+ iterations) finish-step)
61 61
 	      (setf print-flag t)
62 62
 	      (setf stop-flag t))
63
+
63 64
 	    ;; шаг
64 65
 	    (funcall make-step grid newgrid tau)
66
+	    ;; заботимся о выставлении параметров сетки
67
+	    (incf (iteration newgrid))
68
+	    (incf (grid-time newgrid) tau)
65 69
 	    ;; печать
66 70
 	    (when (and (not suppress-export)
67 71
 		       print-flag)
68 72
 	      (export-grid newgrid :print-out-border print-out-border))
69
-	    
73
+
70 74
 	    ;; меняем сетки местами
71 75
 	    (setf grid (clone newgrid))
72 76
 	    ;(update-gradients-conv grid)
... ...
@@ -4,7 +4,6 @@
4 4
   (:use :common-lisp
5 5
 	:alexandria
6 6
 	:gsll)
7
-  
8 7
   (:shadow
9 8
    ;; следующие символы конфликтовали между alexandria и gsll
10 9
    mean variance median standard-deviation factorial
... ...
@@ -142,11 +142,11 @@
142 142
 ;;					  :phi 0.01
143 143
 ;;					  :ug 25
144 144
 ;;					  :ul 400
145
-;;					  :rg 0.001
145
+;;					  :rg 1
146 146
 ;;					  ))
147 147
 ;;
148 148
 ;;(eigen-values-by-cell *eigen-cell* 'A)
149
-;;(eigen-values-by-cell *eigen-cell* 'S)
149
+;;(eigen-values-by-cell *eigen-cell* 'S) => -3e-7
150 150
 ;;(eigen-values-by-cell *eigen-cell* 'Afull)
151 151
 
152 152
 ;; потоковый член в ячейке.
... ...
@@ -357,9 +357,9 @@
357 357
 					 (F3 (/ (* F x4) (* x3 x3)))
358 358
 					 (F4 (/ F x3))
359 359
 					 (TT (/ 1 tau)))
360
-				      (values TT 0 0 0 ;;
360
+				      (values TT 0.0 0.0 0.0 ;;
361 361
 					      F1 (- TT F2) (- F3) F4 ;;
362
-					      0 0 TT 0 ;;
362
+					      0.0 0.0 TT 0.0 ;;
363 363
 					      (- F1) F2 F3 (- TT F4))))
364 364
 				;; собирательная функция (необходима для GSLL)
365 365
 				(eq-system-fdf (x1 x2 x3 x4)
... ...
@@ -370,7 +370,7 @@
370 370
 				      (values v1 v2 v3 v4 j1 j2 j3 j4 j5 j6 j7 j8 j9 j10 j11 j12 j13 j14 j15 j16)))))
371 371
 			 (let-conservative u
372 372
 			   (let ((max-iter *multi-root-iter*)
373
-				 (solver (make-multi-dimensional-root-solver-f 
373
+				 (solver (make-multi-dimensional-root-solver-fdf
374 374
 					  +newton-mfdfsolver+
375 375
 					  (list #'eq-system #'eq-system-df #'eq-system-fdf)
376 376
 					  (grid:make-foreign-array
... ...
@@ -380,9 +380,12 @@
380 380
 					  (or (zerop iter)
381 381
 					      (not (multiroot-test-residual solver *multi-root-est-error*))))
382 382
 				do (iterate solver)
383
-				finally (setf (cell newgrid (index i j))
383
+				finally
384
+				  ;;(format t "solution: ~a~%" (solution solver)
385
+				  (setf (cell newgrid (index i j))
384 386
 					      (make-instance 'conservative-variables
385
-							     :var-vector (solution solver))))))))))
387
+							     :var-vector (vector-double-float->vector
388
+									  (solution solver)))))))))))
386 389
 	   newgrid))
387 390
        (corrector (grid predgrid tau)
388 391
 	 (let ((midgrid (clone grid))
... ...
@@ -408,5 +411,12 @@
408 408
 				    (cvar/ (cvar* tau dF)
409 409
 					   hx))))))
410 410
 	   newgrid)))
411
-    (corrector grid (predictor grid tau) tau)))
411
+    (let ((predgrid (predictor grid tau)))
412
+      (setf newgrid predgrid)
413
+      ;;(corrector grid predgrid tau)
414
+
415
+      (format t "PREDICTOR~%~%")
416
+      (print-grid newgrid :print-out-border t)
417
+      (newline)(newline)
418
+      )))
412 419
 
... ...
@@ -110,7 +110,7 @@
110 110
       )))
111 111
 
112 112
 
113
-(defun test-kt-constant (&key (cells-nx 20) (cells-ny 1) (output-dir #P"~/prog/clfpv/test-KT-constant/"))
113
+(defun test-kt-constant (&key (cells-nx 1) (cells-ny 1) (output-dir #P"~/prog/clfpv/test-KT-constant/"))
114 114
   "Тест на постоянном поле"
115 115
   (ensure-directories-exist output-dir)
116 116
   (let ((*DEFAULT-PATHNAME-DEFAULTS* output-dir))
... ...
@@ -122,22 +122,22 @@
122 122
 			       :bottom-boundary-condition 'zeroe
123 123
 			       :top-boundary-condition 'zeroe)))
124 124
       (init-grid grid
125
-		 :left-border-coord 0
126
-		 :right-border-coord 1
127
-		 :top-border-line (constant-function 1)
128
-		 :bottom-border-line (constant-function 0)
125
+		 :left-border-coord 0.0d0
126
+		 :right-border-coord 1.0d0
127
+		 :top-border-line (constant-function 1.0d0)
128
+		 :bottom-border-line (constant-function 0.0d0)
129 129
 		 :thickening-coef-x 1
130 130
 		 :thickening-coef-y 1)
131 131
       (set-constant-field-2d grid
132 132
 			     (make-instance 'conservative-variables
133
-					    :phi 0.999999
134
-					    :P 100000.0
135
-					    :ug 25.0
136
-					    :ul 25.0))
133
+					    :phi 0.999999d0
134
+					    :P 100000.0d0
135
+					    :ug 25.0d0
136
+					    :ul 25.0d0))
137 137
       (update-boundary-cells grid)
138
-      (march grid :finish-step 1)
139
-      grid
140
-      )))
138
+      (let ((newgrid (march grid :finish-step 1)))
139
+	newgrid
140
+	))))
141 141
 
142 142
 
143 143