/
solver_nr.cpp
7327 lines (6552 loc) · 358 KB
/
solver_nr.cpp
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
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
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
458
459
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
548
549
550
551
552
553
554
555
556
557
558
559
560
561
562
563
564
565
566
567
568
569
570
571
572
573
574
575
576
577
578
579
580
581
582
583
584
585
586
587
588
589
590
591
592
593
594
595
596
597
598
599
600
601
602
603
604
605
606
607
608
609
610
611
612
613
614
615
616
617
618
619
620
621
622
623
624
625
626
627
628
629
630
631
632
633
634
635
636
637
638
639
640
641
642
643
644
645
646
647
648
649
650
651
652
653
654
655
656
657
658
659
660
661
662
663
664
665
666
667
668
669
670
671
672
673
674
675
676
677
678
679
680
681
682
683
684
685
686
687
688
689
690
691
692
693
694
695
696
697
698
699
700
701
702
703
704
705
706
707
708
709
710
711
712
713
714
715
716
717
718
719
720
721
722
723
724
725
726
727
728
729
730
731
732
733
734
735
736
737
738
739
740
741
742
743
744
745
746
747
748
749
750
751
752
753
754
755
756
757
758
759
760
761
762
763
764
765
766
767
768
769
770
771
772
773
774
775
776
777
778
779
780
781
782
783
784
785
786
787
788
789
790
791
792
793
794
795
796
797
798
799
800
801
802
803
804
805
806
807
808
809
810
811
812
813
814
815
816
817
818
819
820
821
822
823
824
825
826
827
828
829
830
831
832
833
834
835
836
837
838
839
840
841
842
843
844
845
846
847
848
849
850
851
852
853
854
855
856
857
858
859
860
861
862
863
864
865
866
867
868
869
870
871
872
873
874
875
876
877
878
879
880
881
882
883
884
885
886
887
888
889
890
891
892
893
894
895
896
897
898
899
900
901
902
903
904
905
906
907
908
909
910
911
912
913
914
915
916
917
918
919
920
921
922
923
924
925
926
927
928
929
930
931
932
933
934
935
936
937
938
939
940
941
942
943
944
945
946
947
948
949
950
951
952
953
954
955
956
957
958
959
960
961
962
963
964
965
966
967
968
969
970
971
972
973
974
975
976
977
978
979
980
981
982
983
984
985
986
987
988
989
990
991
992
993
994
995
996
997
998
999
1000
/* $Id
* Newton-Raphson solver
*/
// ***** Random Useful Code to Note ****
// This code used to be right under the recasting of sol_LU
// Matrix multiplication using superLU-included stuff
// Contents of xvecttest obviously missing
//
//// Try the multiply - X_LU*xvecttest = yvecttest
//// Random attempts - see if we can perform an "Ax+y=y" action
//double *xvecttest, *yvecttest;
//
//// Allocate them
//xvecttest = (double*)gl_malloc(m*sizeof(double));
//yvecttest = (double*)gl_malloc(m*sizeof(double));
//
//sp_dgemv("n",1.0,&A_LU,xvecttest,1,0.0,yvecttest,1);
//// Free, to be thorough
//gl_free(xvecttest);
//gl_free(yvecttest);
// ***** End Random Useful Code to Note ****
/*
***********************************************************************
Generic in-rush notes here
--------------------------
What about P/C relationships for history terms, how handle?
Initialization after returning to service?
***********************************************************************
*/
#include <cmath>
#ifndef GLD_USE_EIGEN
#include "solver_nr.h"
#else
#include "solver_nr_eigen.h"
#endif
/* access to module global variables */
#include "powerflow.h"
#define MT // this enables multithreaded SuperLU
#ifdef MT
#include <slu_mt_ddefs.h> //superLU_MT
#else
#include <slu_ddefs.h> //Sequential superLU (other platforms)
#endif
namespace gld {
template<typename Tv, typename Te>
constexpr Tv pow(Tv value, Te exp) {
switch (exp) {
case 2:
return value * value;
break;
case 3:
return value * value * value;
break;
case 4:
return value * value * value * value;
break;
default:
return std::pow(value, exp);
break;
}
}
}
//SuperLU variable structure
//These are the working variables, but structured for island implementation
typedef struct {
int *perm_c;
int *perm_r;
SuperMatrix A_LU;
SuperMatrix B_LU;
} SUPERLU_NR_vars;
//Initialize the sparse notation
void sparse_init(SPARSE* sm, int nels, int ncols)
{
int indexval;
//Allocate the column pointer GLD heap
sm->cols = (SP_E**)gl_malloc(ncols*sizeof(SP_E*));
//Check it
if (sm->cols == nullptr)
{
GL_THROW("NR: Sparse matrix allocation failed");
/* TROUBLESHOOT
While attempting to allocate space for one of the sparse matrix variables, an error was encountered.
Please try again. If the error persists, please submit your code and a bug report via the ticketing system.
*/
}
else //Zero it, for giggles
{
for (indexval=0; indexval<ncols; indexval++)
{
sm->cols[indexval] = nullptr;
}
}
//Allocate the elements on the GLD heap
sm->llheap = (SP_E*)gl_malloc(nels*sizeof(SP_E));
//Check it
if (sm->llheap == nullptr)
{
GL_THROW("NR: Sparse matrix allocation failed");
//Defined above
}
else //Zero it, for giggles
{
for (indexval=0; indexval<nels; indexval++)
{
sm->llheap[indexval].next = nullptr;
sm->llheap[indexval].row_ind = -1; //Invalid, so it gets upset if we use this (intentional)
sm->llheap[indexval].value = 0.0;
}
}
//Init others
sm->llptr = 0;
sm->ncols = ncols;
}
//Free up/clear the sparse allocations
void sparse_clear(SPARSE* sm)
{
//Clear them up
gl_free(sm->llheap);
gl_free(sm->cols);
//Null them, because I'm paranoid
sm->llheap = nullptr;
sm->cols = nullptr;
//Zero the last one
sm->ncols = 0;
}
void sparse_reset(SPARSE* sm, int ncols)
{
int indexval;
//Set the size
sm->ncols = ncols;
//Do brute force way, for paranoia
for (indexval=0; indexval<ncols; indexval++)
{
sm->cols[indexval] = nullptr;
}
//Set the location pointer
sm->llptr = 0;
}
//Add in new elements to the sparse notation
inline void sparse_add(SPARSE* sm, int row, int col, double value, BUSDATA *bus_values, unsigned int bus_values_count, NR_SOLVER_STRUCT *powerflow_information, int island_number_curr)
{
unsigned int bus_index_val, bus_start_val, bus_end_val;
bool found_proper_bus_val;
SP_E* insertion_point = sm->cols[col];
SP_E* new_list_element = &(sm->llheap[sm->llptr++]);
new_list_element->next = nullptr;
new_list_element->row_ind = row;
new_list_element->value = value;
//if there's a non empty list, traverse to find our rightful insertion point
if(insertion_point != nullptr)
{
if(insertion_point->row_ind > new_list_element->row_ind)
{
//insert the new list element at the first position
new_list_element->next = insertion_point;
sm->cols[col] = new_list_element;
}
else
{
while((insertion_point->next != nullptr) && (insertion_point->next->row_ind < new_list_element->row_ind))
{
insertion_point = insertion_point->next;
}
//Duplicate check -- see how we exited
if (insertion_point->next != nullptr) //We exited because the next element is GEQ to the new element
{
if (insertion_point->next->row_ind == new_list_element->row_ind) //Same entry (by column), so bad
{
//Reset the flag
found_proper_bus_val = false;
//Loop through and see if we can find the bus
for (bus_index_val=0; bus_index_val<bus_values_count; bus_index_val++)
{
//Island check
if (bus_values[bus_index_val].island_number == island_number_curr)
{
//Extract the start/stop indices
bus_start_val = 2*bus_values[bus_index_val].Matrix_Loc;
bus_end_val = bus_start_val + 2*powerflow_information->BA_diag[bus_index_val].size - 1;
//See if we're in this range
if ((new_list_element->row_ind >= bus_start_val) && (new_list_element->row_ind <= bus_end_val))
{
//See if it is actually named -- it should be available
if (bus_values[bus_index_val].name != nullptr)
{
GL_THROW("NR: duplicate admittance entry found - attaches to node %s - check for parallel circuits between common nodes!",bus_values[bus_index_val].name);
/* TROUBLESHOOT
While building up the admittance matrix for the Newton-Raphson solver, a duplicate entry was found.
This is often caused by having multiple lines on the same phases in parallel between two nodes. A reference node
does not have a name. Please name your nodes and and try again. Afterwards, please reconcile this model difference and try again.
*/
}
else
{
GL_THROW("NR: duplicate admittance entry found - no name available - check for parallel circuits between common nodes!");
/* TROUBLESHOOT
While building up the admittance matrix for the Newton-Raphson solver, a duplicate entry was found.
This is often caused by having multiple lines on the same phases in parallel between two nodes. A reference node
does not have a name. Please name your nodes and and try again. Afterwards, please reconcile this model difference and try again.
*/
}//End unnamed bus
}//End found a bus
//Default else -- keep going
}//End part of the island
//Default else -- next bus
}//End of the for loop to find the bus
}//End matches - error
}
else //No next item, so see if our value matches
{
if (insertion_point->row_ind == new_list_element->row_ind) //Same entry (by column), so bad
{
//Reset the flag
found_proper_bus_val = false;
//Loop through and see if we can find the bus
for (bus_index_val=0; bus_index_val<bus_values_count; bus_index_val++)
{
//Island check
if (bus_values[bus_index_val].island_number == island_number_curr)
{
//Extract the start/stop indices
bus_start_val = 2*bus_values[bus_index_val].Matrix_Loc;
bus_end_val = bus_start_val + 2*powerflow_information->BA_diag[bus_index_val].size - 1;
//See if we're in this range
if ((new_list_element->row_ind >= bus_start_val) && (new_list_element->row_ind <= bus_end_val))
{
//See if it is actually named -- it should be available
if (bus_values[bus_index_val].name != nullptr)
{
GL_THROW("NR: duplicate admittance entry found - attaches to node %s - check for parallel circuits between common nodes!",bus_values[bus_index_val].name);
//Defined above
}
else
{
GL_THROW("NR: duplicate admittance entry found - no name available - check for parallel circuits between common nodes!");
//Defined above
}//End unnamed bus
}//End found a bus
//Default else -- keep going
}//End part of the island
//Default else -- next bus
}//End of the for loop to find the bus
}
}
//insert the new list element at the next position
new_list_element->next = insertion_point->next;
insertion_point->next = new_list_element;
}
}
else
sm->cols[col] = new_list_element;
}
void sparse_tonr(SPARSE* sm, NR_SOLVER_VARS *matrices_LUin)
{
//traverse each linked list, which are in order, and copy values into new array
unsigned int rowidx = 0;
unsigned int colidx = 0;
unsigned int i;
SP_E* LL_pointer;
matrices_LUin->cols_LU[0] = 0;
LL_pointer = nullptr;
for(i = 0; i < sm->ncols; i++)
{
LL_pointer = sm->cols[i];
if(LL_pointer != nullptr)
{
matrices_LUin->cols_LU[colidx++] = rowidx;
}
while(LL_pointer != nullptr)
{
matrices_LUin->rows_LU[rowidx] = LL_pointer->row_ind; // row pointers of non zero values
matrices_LUin->a_LU[rowidx] = LL_pointer->value;
++rowidx;
LL_pointer = LL_pointer->next;
}
}
}
/** Newton-Raphson solver
Solves a power flow problem using the Newton-Raphson method
@return n=0 on failure to complete a single iteration,
n>0 to indicate success after n interations, or
n<0 to indicate failure after n iterations
**/
int64 solver_nr(unsigned int bus_count, BUSDATA *bus, unsigned int branch_count, BRANCHDATA *branch, NR_SOLVER_STRUCT *powerflow_values, NRSOLVERMODE powerflow_type , NR_MESHFAULT_IMPEDANCE *mesh_imped_vals, bool *bad_computations)
{
//Index for "islanding operations", when needed
int island_index_val;
//Temporary island looping value
int island_loop_index;
//File pointer for debug outputs
FILE *FPoutVal;
//Saturation mismatch tracking variable
int func_result_val;
//Temporary function variable
FUNCTIONADDR temp_fxn_val;
//Generic status variable
STATUS call_return_status;
//Phase collapser variable
unsigned char phase_worka, phase_workb, phase_workc, phase_workd, phase_worke;
//Temporary calculation variables
double tempIcalcReal, tempIcalcImag;
double tempPbus; //tempPbus store the temporary value of active power load at each bus
double tempQbus; //tempQbus store the temporary value of reactive power load at each bus
//Miscellaneous index variable
unsigned int indexer, tempa, tempb, jindexer, kindexer;
char jindex, kindex;
char temp_index, temp_index_b;
unsigned int temp_index_c;
char mat_temp_index;
//Working matrix for admittance collapsing/determinations
gld::complex tempY[3][3];
//Working matrix for mesh fault impedance storage, prior to "reconstruction"
double temp_z_store[6][6];
//Miscellaneous flag variables
bool Full_Mat_A, Full_Mat_B, proceed_flag, temp_bool_value;
//Deltamode intermediate variables
gld::complex temp_complex_0, temp_complex_1, temp_complex_2, temp_complex_3, temp_complex_4;
gld::complex aval, avalsq;
//Temporary size variable
char temp_size, temp_size_b, temp_size_c;
//Temporary admittance variables
gld::complex Temp_Ad_A[3][3];
gld::complex Temp_Ad_B[3][3];
//DV checking array
gld::complex DVConvCheck[3];
gld::complex currVoltConvCheck[3];
double CurrConvVal;
//Miscellaneous working variable
double work_vals_double_0, work_vals_double_1,work_vals_double_2,work_vals_double_3,work_vals_double_4;
char work_vals_char_0;
bool something_has_been_output;
//SuperLU variables
SuperMatrix L_LU,U_LU;
NCformat *Astore;
DNformat *Bstore;
int nnz, info;
unsigned int m,n;
double *sol_LU;
//Spare notation variable - for output
SP_E *temp_element;
int row, col;
double value;
//Multi-island passes
bool still_iterating_islands;
bool proceed_to_next_island;
int64 return_value_for_solver_NR;
//Multi-island pointer to current superLU variables
SUPERLU_NR_vars *curr_island_superLU_vars;
#ifndef MT
superlu_options_t options; //Additional variables for sequential superLU
SuperLUStat_t stat;
#endif
//Set the global - we're working now, so no more adjustments to island arrays until we're done (except removals)
NR_solver_working = true;
//General "short circuit check" - if there are no islands, just leave
if (NR_islands_detected <= 0)
{
//Make sure the bad computations flag is set
*bad_computations = false;
return 0; //Not really a valid return, but with the false above, should still continue
}
//Ensure bad computations flag is set first
*bad_computations = false;
//Determine special circumstances of SWING bus -- do we want it to truly participate right
if (powerflow_type != PF_NORMAL)
{
if (powerflow_type == PF_DYNCALC) //Parse the list -- anything that is a swing and a generator, deflag it out of principle (and to make it work right)
{
//Set the master swing flag - do for all islands
for (island_loop_index=0; island_loop_index<NR_islands_detected; island_loop_index++)
{
powerflow_values->island_matrix_values[island_loop_index].swing_is_a_swing = false;
}
//Check the buses
for (indexer=0; indexer<bus_count; indexer++)
{
//See if we're a swing-flagged bus
if ((bus[indexer].type > 1) && bus[indexer].swing_functions_enabled)
{
//See if we're "generator ready"
if (*bus[indexer].dynamics_enabled && (bus[indexer].DynCurrent != nullptr))
{
//Deflag us back to "PQ" status
bus[indexer].swing_functions_enabled = false;
//If FPI, flag an admittance update (because now the SWING exists again)
if (NR_solver_algorithm == NRM_FPI)
{
NR_admit_change = true;
}
}
}
//Default else -- normal bus
}//End bus traversion loop
}//Handle running dynamics differently
else //Must be PF_DYNINIT
{
//Loop through the islands
for (island_loop_index=0; island_loop_index<NR_islands_detected; island_loop_index++)
{
//Flag us as true, initially
powerflow_values->island_matrix_values[island_loop_index].swing_is_a_swing = true;
}
}
}//End not normal
else //Must be normal
{
//Loop through all islands
for (island_loop_index=0; island_loop_index<NR_islands_detected; island_loop_index++)
{
powerflow_values->island_matrix_values[island_loop_index].swing_is_a_swing = true; //Flag as a swing, even though this shouldn't do anything
}
}
//Populate aval, if necessary
if (powerflow_type == PF_DYNINIT)
{
//Conversion variables - 1@120-deg
aval = gld::complex(-0.5,(sqrt(3.0)/2.0));
avalsq = aval*aval; //squared value is used a couple places too
}
else //Zero it, just in case something uses it (what would???)
{
aval = 0.0;
avalsq = 0.0;
}
if (matrix_solver_method==MM_EXTERN)
{
for (island_loop_index=0; island_loop_index<NR_islands_detected; island_loop_index++)
{
//Call the initialization routine
powerflow_values->island_matrix_values[island_loop_index].LU_solver_vars = ((void *(*)(void *))(LUSolverFcns.ext_init))(powerflow_values->island_matrix_values[island_loop_index].LU_solver_vars);
//Make sure it worked (allocation check)
if (powerflow_values->island_matrix_values[island_loop_index].LU_solver_vars==nullptr)
{
GL_THROW("External LU matrix solver failed to allocate memory properly!");
/* TROUBLESHOOT
While attempting to allocate memory for the external LU solver, an error occurred.
Please try again. If the error persists, ensure your external LU solver is behaving correctly
and coordinate with their development team as necessary.
*/
}
}//End Island loop
}
else if (matrix_solver_method == MM_SUPERLU) //SuperLU multi-island-related item
{
//Loop through and see if we need to allocate things
for (island_loop_index=0; island_loop_index<NR_islands_detected; island_loop_index++)
{
//See if the superLU variables are populated
if (powerflow_values->island_matrix_values[island_loop_index].LU_solver_vars == nullptr)
{
//Allocate one up
curr_island_superLU_vars = (SUPERLU_NR_vars *)gl_malloc(sizeof(SUPERLU_NR_vars));
//Make sure it worked
if (curr_island_superLU_vars == nullptr)
{
GL_THROW("NR: Failed to allocate memory for one of the necessary matrices");
//Defined elsewhere
}
//Initialize it, for giggles/paranoia
curr_island_superLU_vars->A_LU.Store = nullptr;
//Do the same for the underlying properties, just because
curr_island_superLU_vars->A_LU.Stype = SLU_NC;
curr_island_superLU_vars->A_LU.Dtype = SLU_D;
curr_island_superLU_vars->A_LU.Mtype = SLU_GE;
curr_island_superLU_vars->A_LU.nrow = 0;
curr_island_superLU_vars->A_LU.ncol = 0;
//Repeat for B_LU
curr_island_superLU_vars->B_LU.Store = nullptr;
//Do the same for the underlying properties, just because
curr_island_superLU_vars->B_LU.Stype = SLU_NC;
curr_island_superLU_vars->B_LU.Dtype = SLU_D;
curr_island_superLU_vars->B_LU.Mtype = SLU_GE;
curr_island_superLU_vars->B_LU.nrow = 0;
curr_island_superLU_vars->B_LU.ncol = 0;
//Other arrays
curr_island_superLU_vars->perm_c = nullptr;
curr_island_superLU_vars->perm_r = nullptr;
//Assign it in
powerflow_values->island_matrix_values[island_loop_index].LU_solver_vars = (void *)curr_island_superLU_vars;
}
}
//Renull the random pointer, just in case
curr_island_superLU_vars = nullptr;
}
//If FPI, call the shunt update function
if (NR_solver_algorithm == NRM_FPI)
{
//Bus loop it
for (indexer=0; indexer<bus_count; indexer++)
{
if (bus[indexer].ShuntUpdateFxn != nullptr)
{
//Call the function
call_return_status = ((STATUS (*)(OBJECT *))(*bus[indexer].ShuntUpdateFxn))(bus[indexer].obj);
//Make sure it worked
if (call_return_status == FAILED)
{
GL_THROW("NR: Shunt update failed for device %s",bus[indexer].obj->name ? bus[indexer].obj->name : "Unnamed");
/* TROUBLESHOOT
While attempting to perform the shunt update function call, something failed. Please try again.
If the error persists, please submit your code and a bug report via the ticketing system.
*/
}
//Default else - it worked
}
}//End bus loop
}//End FPI shunt update call
//Call the admittance update code
NR_admittance_update(bus_count,bus,branch_count,branch,powerflow_values,powerflow_type);
//Loop through each island to reset these flags
for (island_loop_index=0; island_loop_index<NR_islands_detected; island_loop_index++)
{
//Reset saturation checks
powerflow_values->island_matrix_values[island_loop_index].SaturationMismatchPresent = false;
//Reset Norton equivalent checks
powerflow_values->island_matrix_values[island_loop_index].NortonCurrentMismatchPresent = false;
//Reset the iteration counters
powerflow_values->island_matrix_values[island_loop_index].iteration_count = 0;
}
//Reset the global flag for iterations on islands -- will be set at the end
still_iterating_islands = true;
//Also reset island index
island_loop_index = 0;
//While it - loop through until we solve the issue
while (still_iterating_islands)
{
//Map the superLU variables each time -- just easier to do it always
if (matrix_solver_method==MM_SUPERLU)
{
//Put the void pointer into a local context
curr_island_superLU_vars = (SUPERLU_NR_vars *)powerflow_values->island_matrix_values[island_loop_index].LU_solver_vars;
}
//Call the load subfunction
compute_load_values(bus_count,bus,powerflow_values,false,island_loop_index);
// Calculate the mismatch of three phase current injection at each bus (deltaI),
//and store the deltaI in terms of real and reactive value in array powerflow_values->island_matrix_values[island_loop_index].current_RHS_NR
if (powerflow_values->island_matrix_values[island_loop_index].current_RHS_NR==nullptr)
{
powerflow_values->island_matrix_values[island_loop_index].current_RHS_NR = (double *)gl_malloc((2*powerflow_values->island_matrix_values[island_loop_index].total_variables) *sizeof(double)); // left_hand side of equation (11)
//Make sure it worked
if (powerflow_values->island_matrix_values[island_loop_index].current_RHS_NR == nullptr)
GL_THROW("NR: Failed to allocate memory for one of the necessary matrices");
//Update the max size
powerflow_values->island_matrix_values[island_loop_index].max_total_variables = powerflow_values->island_matrix_values[island_loop_index].total_variables;
}
else if (powerflow_values->island_matrix_values[island_loop_index].NR_realloc_needed) //Bigger sized (this was checked above)
{
//Decimate the existing value
gl_free(powerflow_values->island_matrix_values[island_loop_index].current_RHS_NR);
//Reallocate it...bigger...faster...stronger!
powerflow_values->island_matrix_values[island_loop_index].current_RHS_NR = (double *)gl_malloc((2*powerflow_values->island_matrix_values[island_loop_index].total_variables) *sizeof(double));
//Make sure it worked
if (powerflow_values->island_matrix_values[island_loop_index].current_RHS_NR == nullptr)
GL_THROW("NR: Failed to allocate memory for one of the necessary matrices");
//Store the updated value
powerflow_values->island_matrix_values[island_loop_index].max_total_variables = powerflow_values->island_matrix_values[island_loop_index].total_variables;
}
//If we're in PF_DYNINIT mode, initialize the balancing convergence
if (powerflow_type == PF_DYNINIT)
{
powerflow_values->island_matrix_values[island_loop_index].swing_converged = true; //init it to yes, fail by exception, not default
}
//Compute the calculated loads (not specified) at each bus
for (indexer=0; indexer<bus_count; indexer++) //for specific bus k
{
//Make sure the bus we're looking at is relevant to this island - otherwise, we can skip it
if (bus[indexer].island_number == island_loop_index)
{
//Mode check
if (NR_solver_algorithm == NRM_FPI)
{
if ((bus[indexer].type > 1) && bus[indexer].swing_functions_enabled)
{
//See if we're a "SWING SWING" or a generator
if (*bus[indexer].dynamics_enabled && (bus[indexer].full_Y != nullptr) && (bus[indexer].DynCurrent != nullptr))
{
//We're not (SWING and SWING-enabled and doing a dynamic init) (not at beginning)
if (!(powerflow_values->island_matrix_values[island_loop_index].swing_is_a_swing && (powerflow_type == PF_DYNINIT)))
{
continue; //Skip us - won't be part of the fixed diagonal
}
}
else
{
continue; //Normal SWING, just skip us
}
}//End SWING and acting as SWING
}//End mode check
//Default TCIM - keep
//Update for generator symmetry - only in static dynamic mode and when SWING is a SWING
if ((powerflow_type == PF_DYNINIT) && powerflow_values->island_matrix_values[island_loop_index].swing_is_a_swing)
{
//Secondary check - see if it is even needed - basically built around 3-phase right now
//Initializes the Norton equivalent item -- normal generators shoudln't need this
if (*bus[indexer].dynamics_enabled && (bus[indexer].full_Y != nullptr) && (bus[indexer].DynCurrent != nullptr))
{
//See if we're three phase or triplex
if ((bus[indexer].phases & 0x07) == 0x07)
{
//Form denominator term of Ii, since it won't change
temp_complex_1 = (~bus[indexer].V[0]) + (~bus[indexer].V[1])*avalsq + (~bus[indexer].V[2])*aval;
//Form up numerator portion that doesn't change (Q and admittance)
//Do in parts, just for readability
//Row 1 of admittance mult
temp_complex_0 = ~bus[indexer].V[0]*(bus[indexer].full_Y[0]*bus[indexer].V[0] + bus[indexer].full_Y[1]*bus[indexer].V[1] + bus[indexer].full_Y[2]*bus[indexer].V[2]);
//Row 2 of admittance
temp_complex_0 += ~bus[indexer].V[1]*(bus[indexer].full_Y[3]*bus[indexer].V[0] + bus[indexer].full_Y[4]*bus[indexer].V[1] + bus[indexer].full_Y[5]*bus[indexer].V[2]);
//Row 3 of admittance
temp_complex_0 += ~bus[indexer].V[2]*(bus[indexer].full_Y[6]*bus[indexer].V[0] + bus[indexer].full_Y[7]*bus[indexer].V[1] + bus[indexer].full_Y[8]*bus[indexer].V[2]);
//Make the conjugate - used for individual phase accumulation later
temp_complex_3 = ~temp_complex_0;
//If we are a SWING, zero our PT portion and QT for accumulation
if (((bus[indexer].type > 1) && bus[indexer].swing_functions_enabled) && (powerflow_values->island_matrix_values[island_loop_index].iteration_count>0))
{
*bus[indexer].PGenTotal = gld::complex(0.0,0.0);
}
}
else if ((bus[indexer].phases & 0x80) == 0x80) //Triplex
{
//Get the "delta voltage" for use here
temp_complex_4 = bus[indexer].V[0] + bus[indexer].V[1];
//Form denominator term of Ii, since it won't change
temp_complex_1 = ~temp_complex_4;
//Form up numerator portion that doesn't change (Q and admittance)
//Should just be a single entry, for "reasons/asssumptions"
temp_complex_0 = ~temp_complex_4*(bus[indexer].full_Y[0]*temp_complex_4);
//Make the conjugate - used for individual phase accumulation later
temp_complex_3 = ~temp_complex_0;
//If we are a SWING, zero our PT portion and QT for accumulation
if (((bus[indexer].type > 1) && bus[indexer].swing_functions_enabled) && (powerflow_values->island_matrix_values[island_loop_index].iteration_count>0))
{
*bus[indexer].PGenTotal = -temp_complex_3; //Gets piecemeal removed from three-phase, just all at once here!
}
}
else if ((bus[indexer].phases & 0x07) == 0x00) //No phases
{
temp_complex_1 = gld::complex(0.0,0.0);
temp_complex_0 = gld::complex(0.0,0.0);
temp_complex_3 = gld::complex(0.0,0.0);
//If we are a SWING, zero our PT portion and QT for accumulation
if (((bus[indexer].type > 1) && bus[indexer].swing_functions_enabled) && (powerflow_values->island_matrix_values[island_loop_index].iteration_count>0))
{
*bus[indexer].PGenTotal = gld::complex(0.0,0.0);
}
}
//Default else - some other permutation, but not really supported (single-phase swing, or partial swing
//numerator done, except PT portion (add in below - SWING bus is different
}
else //Not enabled or not "full-Y-ed" - set to zero
{
temp_complex_0 = 0.0;
temp_complex_1 = 0.0; //Basically a flag to ignore it
}
}//End PF_DYNINIT and SWING is still the same
for (jindex=0; jindex<powerflow_values->BA_diag[indexer].size; jindex++) // rows - for specific phase that exists
{
tempIcalcReal = tempIcalcImag = 0;
if ((bus[indexer].phases & 0x80) == 0x80) //Split phase - triplex bus
{
//Two states of triplex bus - To node of SPCT transformer needs to be different
//First different - Delta-I and diagonal contributions
if ((bus[indexer].phases & 0x20) == 0x20) //We're the To bus
{
//Pre-negated due to the nature of how it's calculated (V1 compared to I1)
tempPbus = bus[indexer].PL[jindex]; //Copy load amounts in
tempQbus = bus[indexer].QL[jindex];
}
else //We're just a normal triplex bus
{
//This one isn't negated (normal operations)
tempPbus = -bus[indexer].PL[jindex]; //Copy load amounts in
tempQbus = -bus[indexer].QL[jindex];
}//end normal triplex bus
if (NR_solver_algorithm == NRM_TCIM)
{
//Get diagonal contributions - only (& always) 2
//Column 1
tempIcalcReal += (powerflow_values->BA_diag[indexer].Y[jindex][0]).Re() * (bus[indexer].V[0]).Re() - (powerflow_values->BA_diag[indexer].Y[jindex][0]).Im() * (bus[indexer].V[0]).Im();// equation (7), the diag elements of bus admittance matrix
tempIcalcImag += (powerflow_values->BA_diag[indexer].Y[jindex][0]).Re() * (bus[indexer].V[0]).Im() + (powerflow_values->BA_diag[indexer].Y[jindex][0]).Im() * (bus[indexer].V[0]).Re();// equation (8), the diag elements of bus admittance matrix
//Column 2
tempIcalcReal += (powerflow_values->BA_diag[indexer].Y[jindex][1]).Re() * (bus[indexer].V[1]).Re() - (powerflow_values->BA_diag[indexer].Y[jindex][1]).Im() * (bus[indexer].V[1]).Im();// equation (7), the diag elements of bus admittance matrix
tempIcalcImag += (powerflow_values->BA_diag[indexer].Y[jindex][1]).Re() * (bus[indexer].V[1]).Im() + (powerflow_values->BA_diag[indexer].Y[jindex][1]).Im() * (bus[indexer].V[1]).Re();// equation (8), the diag elements of bus admittance matrix
}
else //Implies FPI
{
if (powerflow_type == PF_DYNINIT)
{
if ((bus[indexer].type > 1) && bus[indexer].swing_functions_enabled && (bus[indexer].DynCurrent != nullptr)) //We're a generator-type bus
{
//Copy from above, since initialization is same
//Get diagonal contributions - only (& always) 2
//Column 1
tempIcalcReal += (powerflow_values->BA_diag[indexer].Y[jindex][0]).Re() * (bus[indexer].V[0]).Re() - (powerflow_values->BA_diag[indexer].Y[jindex][0]).Im() * (bus[indexer].V[0]).Im();// equation (7), the diag elements of bus admittance matrix
tempIcalcImag += (powerflow_values->BA_diag[indexer].Y[jindex][0]).Re() * (bus[indexer].V[0]).Im() + (powerflow_values->BA_diag[indexer].Y[jindex][0]).Im() * (bus[indexer].V[0]).Re();// equation (8), the diag elements of bus admittance matrix
//Column 2
tempIcalcReal += (powerflow_values->BA_diag[indexer].Y[jindex][1]).Re() * (bus[indexer].V[1]).Re() - (powerflow_values->BA_diag[indexer].Y[jindex][1]).Im() * (bus[indexer].V[1]).Im();// equation (7), the diag elements of bus admittance matrix
tempIcalcImag += (powerflow_values->BA_diag[indexer].Y[jindex][1]).Re() * (bus[indexer].V[1]).Im() + (powerflow_values->BA_diag[indexer].Y[jindex][1]).Im() * (bus[indexer].V[1]).Re();// equation (8), the diag elements of bus admittance matrix
}
//Not a swing-initing-generator, so ignore it
}
//Default else -Not initailization, so no injection if not SWING (and not here)
}
//Now off diagonals
for (kindexer=0; kindexer<(bus[indexer].Link_Table_Size); kindexer++)
{
//Apply proper index to jindexer (easier to implement this way)
jindexer=bus[indexer].Link_Table[kindexer];
if (branch[jindexer].from == indexer) //We're the from bus
{
if ((bus[indexer].phases & 0x20) == 0x20) //SPCT from bus - needs different signage
{
if (NR_solver_algorithm == NRM_TCIM)
{
work_vals_char_0 = jindex*3;
//This situation can only be a normal line (triplex will never be the from for another type)
//Again only, & always 2 columns (just do them explicitly)
//Column 1
tempIcalcReal += ((branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Re() - ((branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Im();// equation (7), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
tempIcalcImag += ((branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Im() + ((branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Re();// equation (8), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
//Column2
tempIcalcReal += ((branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Re() - ((branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Im();// equation (7), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
tempIcalcImag += ((branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Im() + ((branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Re();// equation (8), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
}
else //Implies FPI
{
if (powerflow_type == PF_DYNINIT)
{
if ((bus[indexer].type > 1) && bus[indexer].swing_functions_enabled && (bus[indexer].DynCurrent != nullptr)) //We're a generator-type bus
{
//Copy from above, since initialization is same
work_vals_char_0 = jindex*3;
//This situation can only be a normal line (triplex will never be the from for another type)
//Again only, & always 2 columns (just do them explicitly)
//Column 1
tempIcalcReal += ((branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Re() - ((branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Im();// equation (7), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
tempIcalcImag += ((branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Im() + ((branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Re();// equation (8), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
//Column2
tempIcalcReal += ((branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Re() - ((branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Im();// equation (7), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
tempIcalcImag += ((branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Im() + ((branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Re();// equation (8), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
}
//Not a swing-initing-generator, so ignore it
}
//Not initing, so no care
//SWING bus does something different - injections
if ((bus[branch[jindexer].to].type > 1) && bus[branch[jindexer].to].swing_functions_enabled)
{
work_vals_char_0 = jindex*3;
//This situation can only be a normal line (triplex will never be the from for another type)
//Again only, & always 2 columns (just do them explicitly)
//Column 1
tempIcalcReal -= ((branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Re() - ((branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Im();
tempIcalcImag -= ((branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Im() + ((branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Re();
//Column2
tempIcalcReal -= ((branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Re() - ((branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Im();
tempIcalcImag -= ((branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Im() + ((branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Re();
}
//Default else - normal line, so no injection
}
}//End SPCT To bus - from diagonal contributions
else //Normal line connection to normal triplex
{
if (NR_solver_algorithm == NRM_TCIM)
{
work_vals_char_0 = jindex*3;
//This situation can only be a normal line (triplex will never be the from for another type)
//Again only, & always 2 columns (just do them explicitly)
//Column 1
tempIcalcReal += (-(branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Re() - (-(branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Im();// equation (7), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
tempIcalcImag += (-(branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Im() + (-(branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Re();// equation (8), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
//Column2
tempIcalcReal += (-(branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Re() - (-(branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Im();// equation (7), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
tempIcalcImag += (-(branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Im() + (-(branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Re();// equation (8), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
}
else //Implies FPI
{
if (powerflow_type == PF_DYNINIT)
{
if ((bus[indexer].type > 1) && bus[indexer].swing_functions_enabled && (bus[indexer].DynCurrent != nullptr)) //We're a generator-type bus
{
//Copy from above, since initialization is same
work_vals_char_0 = jindex*3;
//This situation can only be a normal line (triplex will never be the from for another type)
//Again only, & always 2 columns (just do them explicitly)
//Column 1
tempIcalcReal += (-(branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Re() - (-(branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Im();// equation (7), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
tempIcalcImag += (-(branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Im() + (-(branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Re();// equation (8), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
//Column2
tempIcalcReal += (-(branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Re() - (-(branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Im();// equation (7), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
tempIcalcImag += (-(branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Im() + (-(branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Re();// equation (8), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
}
//Not a swing-initing-generator, so ignore it
}
//Not initing, so no care
if ((bus[branch[jindexer].to].type > 1) && bus[branch[jindexer].to].swing_functions_enabled)
{
work_vals_char_0 = jindex*3;
//This situation can only be a normal line (triplex will never be the from for another type)
//Again only, & always 2 columns (just do them explicitly)
//Column 1
tempIcalcReal -= (-(branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Re() - (-(branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Im();
tempIcalcImag -= (-(branch[jindexer].Yfrom[work_vals_char_0])).Re() * (bus[branch[jindexer].to].V[0]).Im() + (-(branch[jindexer].Yfrom[work_vals_char_0])).Im() * (bus[branch[jindexer].to].V[0]).Re();
//Column2
tempIcalcReal -= (-(branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Re() - (-(branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Im();
tempIcalcImag -= (-(branch[jindexer].Yfrom[jindex*3+1])).Re() * (bus[branch[jindexer].to].V[1]).Im() + (-(branch[jindexer].Yfrom[jindex*3+1])).Im() * (bus[branch[jindexer].to].V[1]).Re();
}
//Default else - normal bus, so ignore
}
}//end normal triplex from
}//end from bus
else if (branch[jindexer].to == indexer) //We're the to bus
{
if (branch[jindexer].v_ratio != 1.0) //Transformer
{
//Only a single contributor on the from side - figure out how to get to it
if ((branch[jindexer].phases & 0x01) == 0x01) //C
{
temp_index=2;
}
else if ((branch[jindexer].phases & 0x02) == 0x02) //B
{
temp_index=1;
}
else if ((branch[jindexer].phases & 0x04) == 0x04) //A
{
temp_index=0;
}
else //How'd we get here!?!
{
GL_THROW("NR: A split-phase transformer appears to have an invalid phase");
}
if (NR_solver_algorithm == NRM_TCIM)
{
work_vals_char_0 = jindex*3+temp_index;
//Perform the update, it only happens for one column (nature of the transformer)
tempIcalcReal += (-(branch[jindexer].Yto[work_vals_char_0])).Re() * (bus[branch[jindexer].from].V[temp_index]).Re() - (-(branch[jindexer].Yto[work_vals_char_0])).Im() * (bus[branch[jindexer].from].V[temp_index]).Im();// equation (7), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
tempIcalcImag += (-(branch[jindexer].Yto[work_vals_char_0])).Re() * (bus[branch[jindexer].from].V[temp_index]).Im() + (-(branch[jindexer].Yto[work_vals_char_0])).Im() * (bus[branch[jindexer].from].V[temp_index]).Re();// equation (8), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
}
else //Implies FPI
{
if (powerflow_type == PF_DYNINIT)
{
if ((bus[indexer].type > 1) && bus[indexer].swing_functions_enabled && (bus[indexer].DynCurrent != nullptr)) //We're a generator-type bus
{
//Copy from above, since initialization is same
work_vals_char_0 = jindex*3+temp_index;
//Perform the update, it only happens for one column (nature of the transformer)
tempIcalcReal += (-(branch[jindexer].Yto[work_vals_char_0])).Re() * (bus[branch[jindexer].from].V[temp_index]).Re() - (-(branch[jindexer].Yto[work_vals_char_0])).Im() * (bus[branch[jindexer].from].V[temp_index]).Im();// equation (7), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
tempIcalcImag += (-(branch[jindexer].Yto[work_vals_char_0])).Re() * (bus[branch[jindexer].from].V[temp_index]).Im() + (-(branch[jindexer].Yto[work_vals_char_0])).Im() * (bus[branch[jindexer].from].V[temp_index]).Re();// equation (8), the off_diag elements of bus admittance matrix are equal to negative value of branch admittance
}
//Not a swing-initing-generator, so ignore it
}
//Not initing, so no care
//Bit ridiculous - another "low side triplex is a SWING" scenario
if ((bus[branch[jindexer].from].type > 1) && bus[branch[jindexer].from].swing_functions_enabled)
{
work_vals_char_0 = jindex*3+temp_index;
//Perform the update, it only happens for one column (nature of the transformer)
tempIcalcReal -= (-(branch[jindexer].Yto[work_vals_char_0])).Re() * (bus[branch[jindexer].from].V[temp_index]).Re() - (-(branch[jindexer].Yto[work_vals_char_0])).Im() * (bus[branch[jindexer].from].V[temp_index]).Im();
tempIcalcImag -= (-(branch[jindexer].Yto[work_vals_char_0])).Re() * (bus[branch[jindexer].from].V[temp_index]).Im() + (-(branch[jindexer].Yto[work_vals_char_0])).Im() * (bus[branch[jindexer].from].V[temp_index]).Re();
}
//Default else - other buses need nothing done
}
}//end transformer
else //Must be a normal line then
{
if ((bus[indexer].phases & 0x20) == 0x20) //SPCT from bus - needs different signage
{
if (NR_solver_algorithm == NRM_TCIM)
{
work_vals_char_0 = jindex*3;
//This case should never really exist, but if someone reverses a secondary or is doing meshed secondaries, it might
//Again only, & always 2 columns (just do them explicitly)
//Column 1