📄 spfactor.c
字号:
#else /* Not MODIFIED_MARKOWITZ *//* * QUICK SEARCH OF DIAGONAL FOR PIVOT WITH CONVENTIONAL MARKOWITZ * CRITERION * * Searches the diagonal looking for the best pivot. For a pivot to be * acceptable it must be larger than the pivot RelThreshold times the largest * element in its reduced column. Among the acceptable diagonals, the * one with the smallest MarkowitzProduct is sought. Search terminates * early if a diagonal is found with a MarkowitzProduct of one and its * magnitude is larger than the other elements in its row and column. * Since its MarkowitzProduct is one, there is only one other element in * both its row and column, and, as a condition for early termination, * these elements must be located symmetricly in the matrix. * * >>> Returned: * A pointer to the diagonal element chosen to be pivot. If no diagonal is * acceptable, a NULL is returned. * * >>> Arguments: * Matrix <input> (MatrixPtr) * Pointer to matrix. * Step <input> (int) * Index of the diagonal currently being eliminated. * * >>> Local variables: * ChosenPivot (ElementPtr) * Pointer to the element that has been chosen to be the pivot. * LargestOffDiagonal (RealNumber) * Magnitude of the largest of the off-diagonal terms associated with * a diagonal with MarkowitzProduct equal to one. * Magnitude (RealNumber) * Absolute value of diagonal element. * MinMarkowitzProduct (long) * Smallest Markowitz product found of pivot candidates which are * acceptable. * pDiag (ElementPtr) * Pointer to current diagonal element. * pMarkowitzProduct (long *) * Pointer that points into MarkowitzProduct array. It is used to quickly * access successive Markowitz products. * pOtherInCol (ElementPtr) * When there is only one other element in a column other than the * diagonal, pOtherInCol is used to point to it. Used when Markowitz * product is to determine if off diagonals are placed symmetricly. * pOtherInRow (ElementPtr) * When there is only one other element in a row other than the diagonal, * pOtherInRow is used to point to it. Used when Markowitz product is * to determine if off diagonals are placed symmetricly. */static ElementPtrQuicklySearchDiagonal( Matrix, Step )MatrixPtr Matrix;int Step;{register long MinMarkowitzProduct, *pMarkowitzProduct;register ElementPtr pDiag;int I;ElementPtr ChosenPivot, pOtherInRow, pOtherInCol;RealNumber Magnitude, LargestInCol, LargestOffDiagonal;RealNumber FindBiggestInColExclude();/* Begin `QuicklySearchDiagonal'. */ ChosenPivot = NULL; MinMarkowitzProduct = LARGEST_LONG_INTEGER; pMarkowitzProduct = &(Matrix->MarkowitzProd[Matrix->Size+2]); Matrix->MarkowitzProd[Matrix->Size+1] = Matrix->MarkowitzProd[Step];/* Assure that following while loop will always terminate. */ Matrix->MarkowitzProd[Step-1] = -1;/* * This is tricky. Am using a pointer in the inner while loop to * sequentially step through the MarkowitzProduct array. Search * terminates when the Markowitz product of zero placed at location * Step-1 is found. The row (and column) index on the diagonal is then * calculated by subtracting the pointer to the Markowitz product of * the first diagonal from the pointer to the Markowitz product of the * desired element. The outer for loop is infinite, broken by using * break. * * Search proceeds from the end (high row and column numbers) to the * beginning (low row and column numbers) so that rows and columns with * large Markowitz products will tend to be move to the bottom of the * matrix. However, choosing Diag[Step] is desirable because it would * require no row and column interchanges, so inspect it first by * putting its Markowitz product at the end of the MarkowitzProd * vector. */ for (;;) /* Endless for loop. */ { while (*(--pMarkowitzProduct) >= MinMarkowitzProduct) { /* Just passing through. */ } I = pMarkowitzProduct - Matrix->MarkowitzProd;/* Assure that I is valid; if I < Step, terminate search. */ if (I < Step) break; /* Endless for loop */ if (I > Matrix->Size) I = Step; if ((pDiag = Matrix->Diag[I]) == NULL) continue; /* Endless for loop */ if ((Magnitude = ELEMENT_MAG(pDiag)) <= Matrix->AbsThreshold) continue; /* Endless for loop */ if (*pMarkowitzProduct == 1) {/* Case where only one element exists in row and column other than diagonal. *//* Find off-diagonal elements. */ pOtherInRow = pDiag->NextInRow; pOtherInCol = pDiag->NextInCol; if (pOtherInRow == NULL AND pOtherInCol == NULL) { pOtherInRow = Matrix->FirstInRow[I]; while(pOtherInRow != NULL) { if (pOtherInRow->Col >= Step AND pOtherInRow->Col != I) break; pOtherInRow = pOtherInRow->NextInRow; } pOtherInCol = Matrix->FirstInCol[I]; while(pOtherInCol != NULL) { if (pOtherInCol->Row >= Step AND pOtherInCol->Row != I) break; pOtherInCol = pOtherInCol->NextInCol; } }/* Accept diagonal as pivot if diagonal is larger than off-diagonals and the * off-diagonals are placed symmetricly. */ if (pOtherInRow != NULL AND pOtherInCol != NULL) { if (pOtherInRow->Col == pOtherInCol->Row) { LargestOffDiagonal = MAX(ELEMENT_MAG(pOtherInRow), ELEMENT_MAG(pOtherInCol)); if (Magnitude >= LargestOffDiagonal) {/* Accept pivot, it is unlikely to contribute excess error. */ return pDiag; } } } } MinMarkowitzProduct = *pMarkowitzProduct; ChosenPivot = pDiag; } /* End of endless for loop. */ if (ChosenPivot != NULL) { LargestInCol = FindBiggestInColExclude( Matrix, ChosenPivot, Step ); if( ELEMENT_MAG(ChosenPivot) <= Matrix->RelThreshold * LargestInCol ) ChosenPivot = NULL; } return ChosenPivot;}#endif /* Not MODIFIED_MARKOWITZ *//* * SEARCH DIAGONAL FOR PIVOT WITH MODIFIED MARKOWITZ CRITERION * * Searches the diagonal looking for the best pivot. For a pivot to be * acceptable it must be larger than the pivot RelThreshold times the largest * element in its reduced column. Among the acceptable diagonals, the * one with the smallest MarkowitzProduct is sought. If a tie occurs * between elements of equal MarkowitzProduct, then the element with * the largest ratio between its magnitude and the largest element in its * column is used. The search will be terminated after a given number of * ties have occurred and the best (smallest ratio) of the tied element will * be used as the pivot. The number of ties that will trigger an early * termination is MinMarkowitzProduct * TIES_MULTIPLIER. * * >>> Returned: * A pointer to the diagonal element chosen to be pivot. If no diagonal is * acceptable, a NULL is returned. * * >>> Arguments: * Matrix <input> (MatrixPtr) * Pointer to matrix. * Step <input> (int) * Index of the diagonal currently being eliminated. * * >>> Local variables: * ChosenPivot (ElementPtr) * Pointer to the element that has been chosen to be the pivot. * Size (int) * Local version of size which is placed in a register to increase speed. * Magnitude (RealNumber) * Absolute value of diagonal element. * MinMarkowitzProduct (long) * Smallest Markowitz product found of those pivot candidates which are * acceptable. * NumberOfTies (int) * A count of the number of Markowitz ties that have occurred at current * MarkowitzProduct. * pDiag (ElementPtr) * Pointer to current diagonal element. * pMarkowitzProduct (long*) * Pointer that points into MarkowitzProduct array. It is used to quickly * access successive Markowitz products. * Ratio (RealNumber) * For the current pivot candidate, Ratio is the * Ratio of the largest element in its column to its magnitude. * RatioOfAccepted (RealNumber) * For the best pivot candidate found so far, RatioOfAccepted is the * Ratio of the largest element in its column to its magnitude. */static ElementPtrSearchDiagonal( Matrix, Step )MatrixPtr Matrix;register int Step;{register int J;register long MinMarkowitzProduct, *pMarkowitzProduct;register int I;register ElementPtr pDiag;int NumberOfTies, Size = Matrix->Size;ElementPtr ChosenPivot;RealNumber Magnitude, Ratio, RatioOfAccepted, LargestInCol;RealNumber FindBiggestInColExclude();/* Begin `SearchDiagonal'. */ ChosenPivot = NULL; MinMarkowitzProduct = LARGEST_LONG_INTEGER; pMarkowitzProduct = &(Matrix->MarkowitzProd[Size+2]); Matrix->MarkowitzProd[Size+1] = Matrix->MarkowitzProd[Step];/* Start search of diagonal. */ for (J = Size+1; J > Step; J--) { if (*(--pMarkowitzProduct) > MinMarkowitzProduct) continue; /* for loop */ if (J > Matrix->Size) I = Step; else I = J; if ((pDiag = Matrix->Diag[I]) == NULL) continue; /* for loop */ if ((Magnitude = ELEMENT_MAG(pDiag)) <= Matrix->AbsThreshold) continue; /* for loop *//* Test to see if diagonal's magnitude is acceptable. */ LargestInCol = FindBiggestInColExclude( Matrix, pDiag, Step ); if (Magnitude <= Matrix->RelThreshold * LargestInCol) continue; /* for loop */ if (*pMarkowitzProduct < MinMarkowitzProduct) {/* Notice strict inequality in test. This is a new smallest MarkowitzProduct. */ ChosenPivot = pDiag; MinMarkowitzProduct = *pMarkowitzProduct; RatioOfAccepted = LargestInCol / Magnitude; NumberOfTies = 0; } else {/* This case handles Markowitz ties. */ NumberOfTies++; Ratio = LargestInCol / Magnitude; if (Ratio < RatioOfAccepted) { ChosenPivot = pDiag; RatioOfAccepted = Ratio; } if (NumberOfTies >= MinMarkowitzProduct * TIES_MULTIPLIER) return ChosenPivot; } } /* End of for(Step) */ return ChosenPivot;}#endif /* DIAGONAL_PIVOTING *//* * SEARCH ENTIRE MATRIX FOR BEST PIVOT * * Performs a search over the entire matrix looking for the acceptable * element with the lowest MarkowitzProduct. If there are several that * are tied for the smallest MarkowitzProduct, the tie is broken by using * the ratio of the magnitude of the element being considered to the largest * element in the same column. If no element is acceptable then the largest * element in the reduced submatrix is used as the pivot and the * matrix is declared to be spSMALL_PIVOT. If the largest element is * zero, the matrix is declared to be spSINGULAR. * * >>> Returned: * A pointer to the diagonal element chosen to be pivot. If no element is * found, then NULL is returned and the matrix is spSINGULAR. * * >>> Arguments: * Matrix <input> (MatrixPtr) * Pointer to matrix. * Step <input> (int) * Index of the diagonal currently being eliminated. * * >>> Local variables: * ChosenPivot (ElementPtr) * Pointer to the element that has been chosen to be the pivot. * LargestElementMag (RealNumber) * Magnitude of the largest element yet found in the reduced submatrix. * Size (int) * Local version of Size; placed in a register for speed. * Magnitude (RealNumber) * Absolute value of diagonal element. * MinMarkowitzProduct (long) * Smallest Markowitz product found of pivot candidates which are * acceptable. * NumberOfTies (int) * A count of the number of Markowitz ties that have occurred at current * MarkowitzProduct. * pElement (ElementPtr) * Pointer to current element. * pLargestElement (ElementPtr) * Pointer to the largest element yet found in the reduced submatrix. * Product (long) * Markowitz product for the current row and column. * Ratio (RealNumber) * For the current pivot candidate, Ratio is the * Ratio of the largest element in its column to its magnitude. * RatioOfAccepted (RealNumber) * For the best pivot candidate found so far, RatioOfAccepted is the * Ratio of the largest element in its column to its magnitude. * * >>> Possible errors: * spSINGULAR * spSMALL_PIVOT */static ElementPtrSearchEntireMatrix( Matrix, Step )MatrixPtr Matrix;int Step;{register int I, Size = Matrix->Size;register ElementPtr pElement;int NumberOfTies;long Product, MinMarkowitzProduct;ElementPtr ChosenPivot, pLargestElement;RealNumber Magnitude, LargestElementMag, Ratio, RatioOfAccepted, LargestInCol;RealNumber FindLargestInCol();/* Begin `SearchEntireMatrix'. */ ChosenPivot = NULL; LargestElementMag = 0.0; MinMarkowitzProduct = LARGEST_LONG_INTEGER;/* Start search of matrix on column by column basis. */ for (I = Step; I <= Size; I++) { pElement = Matrix->FirstInCol[I]; while (pElement != NULL AND pElement->Row < Step) pElement = pElement->NextInCol; if((LargestInCol = FindLargestInCol(pElement)) == 0.0) continue; /* for loop */ while (pElement != NULL) {/* Check to see if element is the largest encountered so far. If so, record its magnitude and addr
⌨️ 快捷键说明
复制代码
Ctrl + C
搜索代码
Ctrl + F
全屏模式
F11
切换主题
Ctrl + Shift + D
显示快捷键
?
增大字号
Ctrl + =
减小字号
Ctrl + -