# dfp

a guest
Nov 9th, 2019
90
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
1. /*  Function Rosenbrock
2. F(x,y)=(1-x)^2+100(y-x^2)^2
3. partial derivative x : f = -2+2*x-400*y*x+400*x^3
4. partial derivative y : f = 200*(y-x^2)
5.
6. */
7.
8. #include <conio.h>
9. #include <math.h>
10. #include <iostream>
11.
12.
13.
14. using namespace std;
15.
16.
17.
18. double ** zeroing(double **A) {
19.     for (int i = 0; i < 2; i++) {
20.         for (int j = 0; j < 2; j++)
21.             A[i][j] = 0;
22.     }
23.     return A;
24. }
25.
26. double** product(double **A, double **U) {
27.
28.     double **c = new double*[2];
29.     for (int i = 0; i < 2; i++)
30.         c[i] = new double[2];
31.
32.     for (int i = 0; i < 2; i++) {
33.         for (int j = 0; j < 2; j++) {
34.             c[i][j] = 0;
35.             for (int t = 0; t < 2; t++)
36.                 c[i][j] += A[i][t] * U[t][j];
37.         }
38.     }
39.     return c;
40.
41. }
42.
43. double* product_vector(double **A, double *b) {
44.
45.     double *a = new double[2];
46.
47.     for (int i = 0; i < 2; i++) {
48.         a[i] = 0;
49.         for (int j = 0; j < 2; j++)
50.             a[i] += A[i][j] * b[j];
51.     }
52.     return a;
53. }
54.
55. double  deviation(double **A) {
56.     double sum = 0, div = 0;
57.
58.     for (int i = 0; i < 2; ++i) {
59.         for (int j = 0; j < 2; ++j)
60.             sum += A[i][j] * A[i][j];
61.         div += sum;
62.         sum = 0;
63.     }
64.
65.     return sqrt(div);
66. }
67.
68. int main() {
69.
70.     double k=0, k_max, eps = 0.00001, t, number1, number2;
71.     cout << "enter max approximation: ";
72.     cin >> k_max;
73.     cout << "enter point of min 0<t<1: ";
74.     cin >> t;
75.
76.     double **A = new double*[2]; //identity matrix
77.     for (int i = 0; i < 2; i++)
78.         A[i] = new double[2];
79.
80.     double **B = new double*[2];
81.     for (int i = 0; i < 2; i++)
82.         B[i] = new double[2];
83.
84.     zeroing(B);
85.
86.     double **C = new double*[2];
87.     for (int i = 0; i < 2; i++)
88.         C[i] = new double[2];
89.     zeroing(C);
90.
91.     double **matrix1 = new double*[2];
92.     for (int i = 0; i < 2; i++)
93.         matrix1[i] = new double[2];
94.
95.     zeroing(matrix1);
96.
97.
98.     double **matrix2 = new double*[2];
99.     for (int i = 0; i < 2; i++)
100.         matrix2[i] = new double[2];
101.
102.     zeroing(matrix2);
103.
104.     double *a = new double[2];
105.     for (int i = 0; i < 2; i++)
106.         a[i] = 0;
107.
108.     double *b = new double[2];
109.     for (int i = 0; i < 2; i++)
110.         b[i] = 0;
111.
112.     double *direction = new double[2];
113.     for (int i = 0; i < 2; i++)
114.         direction[i] = 0;
115.
116.     for (int i = 0; i < 2; i++) { // enter identity matrix
117.         for (int j = 0; j < 2; j++) {
118.             if (i == j)
119.                 A[i][j] = 1;
120.             else
121.                 A[i][j] = 0;
122.         }
123.     }
124.     double vector[2] = { 0,0 };
125.     cout << "enter points x and y: ";
126.     cin >> vector[0] >> vector[1];
127.     double gradient[2] = { 0,0 };
128.     double alpha[2] = { 0,0 };
129.     double beta[2] = { 0,0 };
130.     double beta_demo[2] = { 0,0 };
131.
132.
133.     gradient[0] = -2 + 2 * vector[0] - 400 * vector[1] * vector[0] + 400 * pow(vector[0], 3);    // function's gradient evaluation at a point X
134.     gradient[1] = 200 * (vector[1] - pow(vector[0], 2));                                        // function's gradient evaluation at a point Y
136.
137.     do {
138.
139.         number1 = 0, number2 = 0;
140.
142.
143.         /*for (int i = 0; i < 2; i++)
144.             cout << direction[i]<<" ";
145.         cout << endl;*/
146.
147.         direction[0] = -direction[0];
148.         direction[1] = -direction[1];
149.
150.     /*  for (int i = 0; i < 2; i++)
151.             cout << direction[i]<<" ";
152.             cout << endl;*/
153.
154.         for (int i = 0; i < 2; i++)                        // old points
155.             alpha[i] = vector[i];
156.
157.
158.         /*for (int i = 0; i < 2; i++)
159.             cout << alpha[i]<<" ";
160.         cout << endl;
161. */
162.
163.         for (int i = 0; i < 2; i++)                        // new points
164.             vector[i] = vector[i] + t * direction[i];
165.
166.         /*for (int i = 0; i < 2; i++)
167.             cout << vector[i] << " ";
168.         cout << endl;*/
169.
170.         for (int i = 0; i < 2; i++)                        // new points - old points
171.             alpha[i] = vector[i] - alpha[i];
172.
173.     /*  for (int i = 0; i < 2; i++)
174.             cout << alpha[i] << " ";
175.         cout << endl;*/
176.
177.
178.         for (int i = 0; i < 2; i++)                       //old values of grad
180.
181.         /*for (int i = 0; i < 2; i++)
182.             cout << beta[i] << " ";
183.         cout << endl;
184. */
185.
186.         gradient[0] = -2 + 2 * vector[0] - 400 * vector[1] * vector[0] + 400 * pow(vector[0], 3);  //new values of grad
187.         gradient[1] = 200 * (vector[1] - pow(vector[0], 2));
188.
190.
191.         for (int i = 0; i < 2; i++)                         //new grad - old grad
192.             beta[i] = gradient[i] - beta[i];
193.
194.         //for (int i = 0; i < 2; i++)
195.         //  cout << beta[i] << " ";
196.         //cout << endl;
197.
198.
199.
200.         for (int i = 0; i < 2; i++) {                        //begin formula
201.             for (int j = 0; j < 2; j++) {
202.                 matrix1[i][j] = alpha[i] * alpha[j];
203.             }
204.             number1 += alpha[i] * beta[i];
205.         }
206.         /*cout << "vot"<<endl;
207.         for (int i = 0; i < 2; i++) {
208.             for (int j = 0; j < 2; j++)
209.                 cout << matrix1[i][j] << " ";
210.             cout << endl;
211.         }
212.         cout << number1;*/
213.
214.         for (int i = 0; i < 2; i++) {
215.             for (int j = 0; j < 2; j++)
216.                 B[i][j] = matrix1[i][j] / number1;
217.         }
218.
219.
220.         a = product_vector(A, beta);
221.
222.         for (int i = 0; i < 2; i++) {
223.             for (int j = 0; j < 2; j++) {
224.                 matrix2[i][j] = a[i] * beta[j];
225.             }
226.         }
227.
228.         matrix2 = product(matrix2, A);
229.
230.         b = product_vector(A, beta);
231.
232.         for (int i = 0; i < 2; i++) {
233.             number2 += b[i] * beta[i];
234.         }
235.
236.         for (int i = 0; i < 2; i++) {
237.             for (int j = 0; j < 2; j++)
238.                 C[i][j] = matrix2[i][j] / number2;
239.         }
240.
241.         for (int i = 0; i < 2; i++) {
242.             for (int j = 0; j < 2; j++)
243.                 A[i][j]  = A[i][j]+B[i][j]-C[i][j];           //Ak=Ak-1 + B - C;
244.         }
245.         k++;
248.
249.
250.
251.     //cout << sqrt(pow(gradient[0], 2) + pow(gradient[1], 2))<< " " << k<<endl;
252.
253.
254.
255.     cout << vector[0] << " " << vector[1]<<endl;
256.     cout << pow((1 - vector[0]), 2) + 100 * pow((vector[1] - pow(vector[0], 2)), 2);
257.
258.
259.
260.     _getch();
261.     return 0;
262. }
RAW Paste Data