0001 function [u_g varargout]=cielodiagnostic_core_nonlinear(m,params,inputs,analysis)
0002
0003
0004 yts=365*24*3600;
0005
0006
0007 if isempty(inputs),
0008 clear inputs;
0009 end
0010
0011
0012 if isempty(m.y_s),
0013 flag_y_s= 0;
0014 else
0015 flag_y_s= 1;
0016 end
0017
0018
0019 lst=m.lst;
0020
0021
0022 params.kflag=1; params.pflag=1;
0023 params.ktflag=0;
0024
0025
0026 count=1;
0027 converged=0;
0028 soln(count).u_g=[];
0029 soln(count).u_f=[];
0030
0031 if params.debug,
0032 disp(sprintf('\n%s',[' starting direct shooting method']));
0033 end
0034 while(~converged),
0035
0036
0037
0038 inputs.velocity=full(soln(count).u_g);
0039
0040 [K_gg_nopenalty , p_g_nopenalty , dK_gg_nopenalty]=Emg(m.bgpdt,m.bgpdtb, m.est,lst,m.ept,m.mpt,m.geom3,params,inputs,analysis);
0041 [K_gg , p_g]=PenaltyEmg(K_gg_nopenalty,p_g_nopenalty,m.bgpdt,m.bgpdtb, m.est,lst,m.ept,m.mpt,m.geom3,params,inputs,analysis);
0042
0043
0044 [K_ff, K_fs] = Reducematrixfromgtof( K_gg, m.G_mn, flag_y_s ,m.uset);
0045
0046
0047 [p_f] = Reducerightside( p_g, m.G_mn, K_fs, m.y_s, flag_y_s ,m.uset);
0048
0049
0050 count=count+1;
0051
0052
0053 [soln(count).u_f]=Solver(K_ff,p_f,'general');
0054
0055
0056 [soln(count).u_g]= Mergesolvec( soln(count).u_f, m.G_mn, m.y_s, m.uset );
0057
0058 if (count>2),
0059 soln(count-1).u_f=NaN;
0060 end
0061
0062
0063 inputs.velocity=soln(count).u_g;
0064
0065 [lst,constraints_converged,num_unstable_constraints] =PenaltyConstraints( m.bgpdt,m.bgpdtb, m.est,lst,m.mpt,inputs,analysis);
0066
0067
0068 if(count>=3),
0069 dug=soln(count).u_g-soln(count-1).u_g; ndu=norm(dug,2); nu=norm(soln(count-1).u_g,2);
0070 if (ndu/nu<=params.eps_rel),
0071 if constraints_converged,
0072 converged=1;
0073 end
0074 end
0075
0076 if params.debug,
0077 disp(sprintf('%s %g %s %g',' convergence criterion: norm(du)/norm(u)=',ndu/nu,' > ',params.eps_rel));
0078 disp(sprintf('%s %i',' number of unstable constraints=',num_unstable_constraints));
0079 end
0080 end
0081
0082 if (count>2),
0083 soln(count-1).u_g=NaN;
0084 end
0085 end
0086
0087
0088 [soln(count).u_f]=NaN;
0089
0090 u_g=soln(count).u_g;
0091
0092
0093 nout=max(nargout,1)-1;
0094 if nout==2,
0095
0096 inputs.velocity=full(soln(count).u_g);
0097 params.kflag=1; params.pflag=0; params.ktflag=0;
0098 [soln(count).K_gg , soln(count).p_g , soln(count).dK_gg]=Emg(m.bgpdt,m.bgpdtb, m.est,lst,m.ept,m.mpt,m.geom3,params,inputs,analysis);
0099 [soln(count).K_ff, soln(count).K_fs] = Reducematrixfromgtof( soln(count).K_gg, m.G_mn, flag_y_s, m.uset );
0100 varargout(1)={soln(count).K_ff};
0101 varargout(2)={soln(count).K_fs};
0102 end
0103
0104 end