Overview
Comment:Modified step size adjustment
Downloads: Tarball | ZIP archive
Timelines: family | ancestors | descendants | both | origin/master | trunk
Files: files | file ages | folders
SHA3-256: 4c14f1758d365dae57cef3a4465e3f8b56b09aed76e7e75e3c8bc5128e4661f2
User & Date: gawthrop@users.sourceforge.net on 1999-11-11 21:48:46.000
Other Links: branch diff | manifest | tags
Context
1999-11-12
06:56:18
*** empty log message *** check-in: cd23f4d0ce user: gawthrop@users.sourceforge.net tags: origin/master, trunk
1999-11-11
21:48:46
Modified step size adjustment check-in: 4c14f1758d user: gawthrop@users.sourceforge.net tags: origin/master, trunk
1999-11-10
00:47:08
Replaced ifs by a table of cr/arg information check-in: b0b149b318 user: gawthrop@users.sourceforge.net tags: origin/master, trunk
Changes
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43


44
45
46
47
48
49
50
51
52
53
54

55
56
















57
58
59



60

61
62
63


64
65
66

67
68
69

70
71
72
  endif
  
  if nargin<6
    weight = ones(size(y_s));
  endif

  if nargin<7
    criterion = 1e-5;
  endif
  
  if nargin <8
    max_iterations = 10;
  endif

  if nargin<9
    alpha = 1.0;
  endif
  
  if (!strcmp(method,"time"))&&(!strcmp(method,"freq"))
    error("method must be either time or freq")
  endif
  


  N_theta = length(free);
  Weight = weight*ones(1,N_theta); # Sensitivity weight
  e_last = 1e20;
  error=1e10;
  theta = theta_0;
  Theta = [];
  Error = [];
  Y = [];
  iterations = -1;
  while (abs(e_last-error)>criterion)&&(iterations<max_iterations)
    iterations = iterations + 1;

    e_last = error;
    eval(sprintf("[t,y,y_theta] = mtt_s%s(system_name,theta,free);",method)); # Simulate system
















    Theta = [Theta theta];	# Save parameters
    Y = [Y y];			# Save output
    E = weight.*(y - y_s);	# Weighted error



    y_theta = Weight.*y_theta;	# Weighted sensitivity

    error = (E'*E);		# Sum the error
    Error = [Error error];
    ##    theta(free) = theta(free) - alpha*(real(y_theta'*y_theta)\real(y_theta'*E));


    tol = 1e-4;
    JJ = real(y_theta'*y_theta);
    ## sigma = svd(JJ)

    theta(free) = theta(free) - alpha*( pinv(JJ,tol)*real(y_theta'*E) );
  endwhile


endfunction









|



|










>
>
|
<
|






|

>
|

>
>
>
>
>
>
>
>
>
>
>
>
>
>
>
>


|
>
>
>
|
>
|
|
|
>
>
|
|
|
>
|
|

>



22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46

47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
  endif
  
  if nargin<6
    weight = ones(size(y_s));
  endif

  if nargin<7
    criterion = 1e-7;
  endif
  
  if nargin <8
    max_iterations = 25;
  endif

  if nargin<9
    alpha = 1.0;
  endif
  
  if (!strcmp(method,"time"))&&(!strcmp(method,"freq"))
    error("method must be either time or freq")
  endif
  
  [n_data,n_y] = size(y_s);

  n_th = length(free);

  error_old = 1e20;
  error=1e10;
  theta = theta_0;
  Theta = [];
  Error = [];
  Y = [];
  iterations = -1;
  while (abs(error_old-error)>criterion)&&(iterations<max_iterations)
    iterations = iterations + 1;
    error_old_old = error_old;
    error_old = error;
    eval(sprintf("[t,y,y_theta] = mtt_s%s(system_name,theta,free);",method)); # Simulate system
    error = 0; 
    J = zeros(n_th,1);
    JJ = zeros(n_th,n_th);
    for i = 1:n_y
      E = weight(:,i).*(y(:,i) - y_s(:,i)); # Weighted error
      error = error + (E'*E);	# Sum the error

      Weight = weight(:,i)*ones(1,n_th); # Sensitivity weight
      y_theta_w = Weight.*y_theta(:,i:n_y:n_y*n_th); # Weighted sensitivity
      J  = J + real(y_theta_w'*E); # Jacobian
      JJ = JJ + real(y_theta_w'*y_theta_w); # Newton Euler approx Hessian
    endfor

    error = error
    ## Diagnostics
    Error = [Error error];	# Save error
    Theta = [Theta theta];	# Save parameters
    Y = [Y y];			# Save output

    ## Update the estimate if we are not done yet.
    if (abs(error_old-error)>criterion)&&(iterations<max_iterations)
      if error>error_old+criterion # Halve step size and try again
        factor = 10;
	disp(sprintf("step/%i",factor));
	error = error_old;	# Go back to previous error value
	error_old = inf;	# Don't let it think its converged
	theta(free) = theta(free)  + step; # Reverse step
	step = step/factor;	# new stepsize
      else			# Recompute step size
	tol = 1e-5;
	step =  pinv(JJ,tol)*J;	# Step size
	#step =  pinv(JJ)*J;	# Step size (built in tol)
      endif
      theta(free) = theta(free) - step; # Increment parameters
    endif

   endwhile
endfunction



MTT: Model Transformation Tools
GitHub | SourceHut | Sourceforge | Fossil RSS ]