(self, rhs : BaseVector, sol : BaseVector)
| 572 | super().__init__(*args, **kwargs) |
| 573 | |
| 574 | def _SolveImpl(self, rhs : BaseVector, sol : BaseVector): |
| 575 | # following implementation from scipy |
| 576 | mat,pre,tol = self.mat, self.pre, self.tol |
| 577 | |
| 578 | r = rhs.CreateVector() |
| 579 | u = rhs.CreateVector() |
| 580 | v = rhs.CreateVector() |
| 581 | w = rhs.CreateVector() |
| 582 | uhat = rhs.CreateVector() |
| 583 | uNext = rhs.CreateVector() |
| 584 | rstar = rhs.CreateVector() |
| 585 | d = rhs.CreateVector() |
| 586 | x = rhs.CreateVector() |
| 587 | z = rhs.CreateVector() |
| 588 | tmp = rhs.CreateVector() |
| 589 | d[:] = 0 |
| 590 | |
| 591 | if Norm(rhs)==0: |
| 592 | sol[:] = 0 |
| 593 | return |
| 594 | |
| 595 | x.data = sol |
| 596 | r.data = rhs - mat*sol |
| 597 | |
| 598 | u.data = r |
| 599 | w.data = r |
| 600 | rstar.data = r |
| 601 | |
| 602 | v.data = pre@mat * r |
| 603 | uhat.data = v |
| 604 | |
| 605 | theta = eta = 0 |
| 606 | |
| 607 | rho = InnerProduct(rstar, r) |
| 608 | rhoLast = rho |
| 609 | r0norm = sqrt(rho) |
| 610 | |
| 611 | tau = r0norm |
| 612 | if r0norm == 0: |
| 613 | return |
| 614 | |
| 615 | for iter in range(0,self.maxiter): |
| 616 | even = iter%2 == 0 |
| 617 | if (even): |
| 618 | vtrstar = InnerProduct(rstar, v) |
| 619 | if vtrstar==0: |
| 620 | return |
| 621 | |
| 622 | alpha = rho/vtrstar |
| 623 | uNext.data = u-alpha*v |
| 624 | |
| 625 | w -= alpha*uhat |
| 626 | d *= theta**2/alpha*eta |
| 627 | d += u |
| 628 | theta = Norm(w)/tau |
| 629 | c = sqrt(1/(1+theta**2)) |
| 630 | tau *= theta*c |
| 631 |
nothing calls this directly
no test coverage detected