1 | function femmodel=solver_couplingstokes_nonlinear(femmodel,conserve_loads)
|
---|
2 | %SOLVER_COUPLINGSTOKES_NONLINEAR - core solver of coupling run
|
---|
3 | %
|
---|
4 | % Usage:
|
---|
5 | % [femmodel]=solver_couplingstokes_nonlinear(femmodel,conserve_loads)
|
---|
6 |
|
---|
7 | %initialize solution vector
|
---|
8 | converged=0; count=1;
|
---|
9 |
|
---|
10 | %First get ug=ug_horiz+ug_vert
|
---|
11 | femmodel=SetCurrentConfiguration(femmodel,DiagnosticHorizAnalysisEnum);
|
---|
12 | ug_horiz=GetSolutionFromInputs(femmodel.elements, femmodel.nodes, femmodel.vertices, femmodel.loads, femmodel.materials, femmodel.parameters);
|
---|
13 | uf_horiz=Reducevectorgtof( ug_horiz, femmodel.nodes,femmodel.parameters);
|
---|
14 |
|
---|
15 | while(~converged),
|
---|
16 |
|
---|
17 | %First compute the horizontal velocity
|
---|
18 | femmodel=SetCurrentConfiguration(femmodel,DiagnosticHorizAnalysisEnum);
|
---|
19 | configuration_type=femmodel.parameters.ConfigurationType;
|
---|
20 |
|
---|
21 | %Update the solution to make sure that vx and vxold are similar
|
---|
22 | [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,ug_horiz);
|
---|
23 |
|
---|
24 | %save pointer to old velocity
|
---|
25 | old_ug_horiz=ug_horiz;
|
---|
26 | old_uf_horiz=uf_horiz;
|
---|
27 |
|
---|
28 | [K_ff_horiz,K_fs_horiz,p_f_horiz,d_f_horiz,kmax_horiz]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
|
---|
29 | ys=CreateNodalConstraints(femmodel.nodes,configuration_type);
|
---|
30 | p_f_horiz = Reduceload( p_f_horiz, K_fs_horiz, ys);
|
---|
31 |
|
---|
32 | uf_horiz=Solver(K_ff_horiz,p_f_horiz,old_uf_horiz,d_f_horiz,femmodel.parameters);
|
---|
33 | ug_horiz= Mergesolutionfromftog( uf_horiz, ys, femmodel.nodes,femmodel.parameters);
|
---|
34 |
|
---|
35 | [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,ug_horiz);
|
---|
36 |
|
---|
37 | %Figure out if convergence have been reached
|
---|
38 | converged=convergence(K_ff_horiz,p_f_horiz,uf_horiz,old_uf_horiz,femmodel.parameters);
|
---|
39 |
|
---|
40 | %Then compute vertical velocity
|
---|
41 | femmodel=SetCurrentConfiguration(femmodel,DiagnosticVertAnalysisEnum);
|
---|
42 | configuration_type=femmodel.parameters.ConfigurationType;
|
---|
43 |
|
---|
44 | [K_ff_vert,K_fs_vert,p_f_vert,d_f_vert,kmax_vert]=SystemMatrices(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters);
|
---|
45 | ys=CreateNodalConstraints(femmodel.nodes,configuration_type);
|
---|
46 | p_f_vert = Reduceload( p_f_vert, K_fs_vert, ys);
|
---|
47 |
|
---|
48 | uf_vert=Solver(K_ff_vert,p_f_vert,[],d_f_vert,femmodel.parameters);
|
---|
49 | ug_vert= Mergesolutionfromftog( uf_vert, ys, femmodel.nodes,femmodel.parameters);
|
---|
50 |
|
---|
51 | [femmodel.elements,femmodel.materials]=InputUpdateFromSolution(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,ug_vert);
|
---|
52 |
|
---|
53 | end
|
---|
54 | end
|
---|