home *** CD-ROM | disk | FTP | other *** search
Unknown | 1994-01-13 | 7.6 KB |
open in:
MacOS 8.1
|
Win98
|
DOS
view JSON data
|
view as text
This file was not able to be converted.
This format is not currently supported by dexvert.
Confidence | Program | Detection | Match Type | Support
|
---|
100%
| file
| data
| default
|
|
hex view+--------+-------------------------+-------------------------+--------+--------+
|00000000| 65 72 66 61 63 65 20 74 | 6f 20 69 74 65 72 5f 63 |erface t|o iter_c|
|00000010| 67 6e 65 28 29 20 77 68 | 69 63 68 20 0a 20 20 20 |gne() wh|ich . |
|00000020| 75 73 65 73 20 73 70 61 | 72 73 65 20 6d 61 74 72 |uses spa|rse matr|
|00000030| 69 78 20 64 61 74 61 20 | 73 74 72 75 63 74 75 72 |ix data |structur|
|00000040| 65 73 0a 20 20 20 2d 2d | 20 61 73 73 75 6d 65 73 |es. --| assumes|
|00000050| 20 74 68 61 74 20 42 20 | 63 6f 6e 74 61 69 6e 73 | that B |contains|
|00000060| 20 61 6e 20 61 63 74 75 | 61 6c 20 70 72 65 63 6f | an actu|al preco|
|00000070| 6e 64 69 74 69 6f 6e 65 | 72 20 28 6f 72 20 4e 55 |nditione|r (or NU|
|00000080| 4c 4c 29 0a 20 20 20 75 | 73 65 20 61 6c 77 61 79 |LL). u|se alway|
|00000090| 73 20 61 73 20 66 6f 6c | 6c 6f 77 73 3a 0a 20 20 |s as fol|lows:. |
|000000a0| 20 20 20 20 78 20 3d 20 | 69 74 65 72 5f 73 70 63 | x = |iter_spc|
|000000b0| 67 6e 65 28 41 2c 42 2c | 62 2c 65 70 73 2c 78 2c |gne(A,B,|b,eps,x,|
|000000c0| 6c 69 6d 69 74 2c 73 74 | 65 70 73 29 3b 0a 20 20 |limit,st|eps);. |
|000000d0| 20 6f 72 20 0a 20 20 20 | 20 20 20 78 20 3d 20 69 | or . | x = i|
|000000e0| 74 65 72 5f 73 70 63 67 | 6e 65 28 41 2c 42 2c 62 |ter_spcg|ne(A,B,b|
|000000f0| 2c 65 70 73 2c 56 4e 55 | 4c 4c 2c 6c 69 6d 69 74 |,eps,VNU|LL,limit|
|00000100| 2c 73 74 65 70 73 29 3b | 0a 20 20 20 49 6e 20 74 |,steps);|. In t|
|00000110| 68 65 20 73 65 63 6f 6e | 64 20 63 61 73 65 20 74 |he secon|d case t|
|00000120| 68 65 20 73 6f 6c 75 74 | 69 6f 6e 20 76 65 63 74 |he solut|ion vect|
|00000130| 6f 72 20 69 73 20 63 72 | 65 61 74 65 64 2e 0a 2a |or is cr|eated..*|
|00000140| 2f 0a 56 45 43 20 20 2a | 69 74 65 72 5f 73 70 63 |/.VEC *|iter_spc|
|00000150| 67 6e 65 28 41 2c 42 2c | 62 2c 65 70 73 2c 78 2c |gne(A,B,|b,eps,x,|
|00000160| 6c 69 6d 69 74 2c 73 74 | 65 70 73 29 0a 53 50 4d |limit,st|eps).SPM|
|00000170| 41 54 09 2a 41 2c 20 2a | 42 3b 0a 56 45 43 09 2a |AT.*A, *|B;.VEC.*|
|00000180| 62 2c 20 2a 78 3b 0a 64 | 6f 75 62 6c 65 09 65 70 |b, *x;.d|ouble.ep|
|00000190| 73 3b 0a 69 6e 74 20 2a | 73 74 65 70 73 2c 20 6c |s;.int *|steps, l|
|000001a0| 69 6d 69 74 3b 0a 7b 09 | 0a 20 20 20 49 54 45 52 |imit;.{.|. ITER|
|000001b0| 20 2a 69 70 3b 0a 20 20 | 20 0a 20 20 20 69 70 20 | *ip;. | . ip |
|000001c0| 3d 20 69 74 65 72 5f 67 | 65 74 28 30 2c 30 29 3b |= iter_g|et(0,0);|
|000001d0| 0a 20 20 20 69 70 2d 3e | 41 78 20 3d 20 28 46 75 |. ip->|Ax = (Fu|
|000001e0| 6e 5f 41 78 29 20 73 70 | 5f 6d 76 5f 6d 6c 74 3b |n_Ax) sp|_mv_mlt;|
|000001f0| 0a 20 20 20 69 70 2d 3e | 41 5f 70 61 72 20 3d 20 |. ip->|A_par = |
|00000200| 28 76 6f 69 64 20 2a 29 | 41 3b 0a 20 20 20 69 70 |(void *)|A;. ip|
|00000210| 2d 3e 41 54 78 20 3d 20 | 28 46 75 6e 5f 41 78 29 |->ATx = |(Fun_Ax)|
|00000220| 20 73 70 5f 76 6d 5f 6d | 6c 74 3b 0a 20 20 20 69 | sp_vm_m|lt;. i|
|00000230| 70 2d 3e 41 54 5f 70 61 | 72 20 3d 20 28 76 6f 69 |p->AT_pa|r = (voi|
|00000240| 64 20 2a 29 41 3b 0a 20 | 20 20 69 66 20 28 42 29 |d *)A;. | if (B)|
|00000250| 20 7b 0a 20 20 20 20 20 | 20 69 70 2d 3e 42 78 20 | {. | ip->Bx |
|00000260| 3d 20 28 46 75 6e 5f 41 | 78 29 20 73 70 5f 6d 76 |= (Fun_A|x) sp_mv|
|00000270| 5f 6d 6c 74 3b 0a 20 20 | 20 20 20 20 69 70 2d 3e |_mlt;. | ip->|
|00000280| 42 5f 70 61 72 20 3d 20 | 28 76 6f 69 64 20 2a 29 |B_par = |(void *)|
|00000290| 42 3b 0a 20 20 20 7d 0a | 20 20 20 65 6c 73 65 20 |B;. }.| else |
|000002a0| 7b 0a 20 20 20 20 20 20 | 69 70 2d 3e 42 78 20 3d |{. |ip->Bx =|
|000002b0| 20 28 46 75 6e 5f 41 78 | 29 20 4e 55 4c 4c 3b 0a | (Fun_Ax|) NULL;.|
|000002c0| 20 20 20 20 20 20 69 70 | 2d 3e 42 5f 70 61 72 20 | ip|->B_par |
|000002d0| 3d 20 4e 55 4c 4c 3b 0a | 20 20 20 7d 0a 20 20 20 |= NULL;.| }. |
|000002e0| 69 70 2d 3e 69 6e 66 6f | 20 3d 20 28 46 75 6e 5f |ip->info| = (Fun_|
|000002f0| 69 6e 66 6f 29 20 4e 55 | 4c 4c 3b 0a 20 20 20 69 |info) NU|LL;. i|
|00000300| 70 2d 3e 62 20 3d 20 62 | 3b 0a 20 20 20 69 70 2d |p->b = b|;. ip-|
|00000310| 3e 65 70 73 20 3d 20 65 | 70 73 3b 0a 20 20 20 69 |>eps = e|ps;. i|
|00000320| 70 2d 3e 6c 69 6d 69 74 | 20 3d 20 6c 69 6d 69 74 |p->limit| = limit|
|00000330| 3b 0a 20 20 20 69 70 2d | 3e 78 20 3d 20 78 3b 0a |;. ip-|>x = x;.|
|00000340| 20 20 20 69 74 65 72 5f | 63 67 6e 65 28 69 70 29 | iter_|cgne(ip)|
|00000350| 3b 0a 20 20 20 78 20 3d | 20 69 70 2d 3e 78 3b 0a |;. x =| ip->x;.|
|00000360| 20 20 20 69 66 20 28 73 | 74 65 70 73 29 20 2a 73 | if (s|teps) *s|
|00000370| 74 65 70 73 20 3d 20 69 | 70 2d 3e 73 74 65 70 73 |teps = i|p->steps|
|00000380| 3b 0a 20 20 20 69 70 2d | 3e 73 68 61 72 65 64 5f |;. ip-|>shared_|
|00000390| 78 20 3d 20 69 70 2d 3e | 73 68 61 72 65 64 5f 62 |x = ip->|shared_b|
|000003a0| 20 3d 20 54 52 55 45 3b | 0a 20 20 20 69 74 65 72 | = TRUE;|. iter|
|000003b0| 5f 66 72 65 65 28 69 70 | 29 3b 20 20 20 2f 2a 20 |_free(ip|); /* |
|000003c0| 72 65 6c 65 61 73 65 20 | 6f 6e 6c 79 20 49 54 45 |release |only ITE|
|000003d0| 52 20 73 74 72 75 63 74 | 75 72 65 20 2a 2f 0a 20 |R struct|ure */. |
|000003e0| 20 20 72 65 74 75 72 6e | 20 78 3b 09 09 0a 7d 0a | return| x;...}.|
|000003f0| 0a 0a 0a 0a 46 69 6c 65 | 44 61 74 61 82 69 74 65 |....File|Data.ite|
|00000400| 72 73 79 6d 00 00 00 00 | 02 36 00 00 45 ff ff ff |rsym....|.6..E...|
|00000410| 38 53 40 17 03 00 56 b4 | 0a 2f 2a 2a 2a 2a 2a 2a |8S@...V.|./******|
|00000420| 2a 2a 2a 2a 2a 2a 2a 2a | 2a 2a 2a 2a 2a 2a 2a 2a |********|********|
|00000430| 2a 2a 2a 2a 2a 2a 2a 2a | 2a 2a 2a 2a 2a 2a 2a 2a |********|********|
|00000440| 2a 2a 2a 2a 2a 2a 2a 2a | 2a 2a 2a 2a 2a 2a 2a 2a |********|********|
|00000450| 2a 2a 2a 2a 2a 2a 2a 2a | 2a 2a 2a 2a 2a 2a 2a 2a |********|********|
|00000460| 2a 2a 2a 2a 0a 2a 2a 0a | 2a 2a 20 43 6f 70 79 72 |****.**.|** Copyr|
|00000470| 69 67 68 74 20 28 43 29 | 20 31 39 39 33 20 44 61 |ight (C)| 1993 Da|
|00000480| 76 69 64 20 45 2e 20 53 | 74 65 77 61 72 64 20 26 |vid E. S|teward &|
|00000490| 20 5a 62 69 67 6e 69 65 | 77 20 4c 65 79 6b 2c 20 | Zbignie|w Leyk, |
|000004a0| 61 6c 6c 20 72 69 67 68 | 74 73 20 72 65 73 65 72 |all righ|ts reser|
|000004b0| 76 65 64 2e 0a 2a 2a 0a | 2a 2a 09 09 09 20 20 20 |ved..**.|**... |
|000004c0| 20 20 4d 65 73 63 68 61 | 63 68 20 4c 69 62 72 61 | Mescha|ch Libra|
|000004d0| 72 79 0a 2a 2a 20 0a 2a | 2a 20 54 68 69 73 20 4d |ry.** .*|* This M|
|000004e0| 65 73 63 68 61 63 68 20 | 4c 69 62 72 61 72 79 20 |eschach |Library |
|000004f0| 69 73 20 70 72 6f 76 69 | 64 65 64 20 22 61 73 20 |is provi|ded "as |
|00000500| 69 73 22 20 77 69 74 68 | 6f 75 74 20 61 6e 79 20 |is" with|out any |
|00000510| 65 78 70 72 65 73 73 20 | 0a 2a 2a 20 6f 72 20 69 |express |.** or i|
|00000520| 6d 70 6c 69 65 64 20 77 | 61 72 72 61 6e 74 79 20 |mplied w|arranty |
|00000530| 6f 66 20 61 6e 79 20 6b | 69 6e 64 20 77 69 74 68 |of any k|ind with|
|00000540| 20 72 65 73 70 65 63 74 | 20 74 6f 20 74 68 69 73 | respect| to this|
|00000550| 20 73 6f 66 74 77 61 72 | 65 2e 20 0a 2a 2a 20 49 | softwar|e. .** I|
|00000560| 6e 20 70 61 72 74 69 63 | 75 6c 61 72 20 74 68 65 |n partic|ular the|
|00000570| 20 61 75 74 68 6f 72 73 | 20 73 68 61 6c 6c 20 6e | authors| shall n|
|00000580| 6f 74 20 62 65 20 6c 69 | 61 62 6c 65 20 66 6f 72 |ot be li|able for|
|00000590| 20 61 6e 79 20 64 69 72 | 65 63 74 2c 20 0a 2a 2a | any dir|ect, .**|
|000005a0| 20 69 6e 64 69 72 65 63 | 74 2c 20 73 70 65 63 69 | indirec|t, speci|
|000005b0| 61 6c 2c 20 69 6e 63 69 | 64 65 6e 74 61 6c 20 6f |al, inci|dental o|
|000005c0| 72 20 63 6f 6e 73 65 71 | 75 65 6e 74 69 61 6c 20 |r conseq|uential |
|000005d0| 64 61 6d 61 67 65 73 20 | 61 72 69 73 69 6e 67 20 |damages |arising |
|000005e0| 0a 2a 2a 20 69 6e 20 61 | 6e 79 20 77 61 79 20 66 |.** in a|ny way f|
|000005f0| 72 6f 6d 20 75 73 65 20 | 6f 66 20 74 68 65 20 73 |rom use |of the s|
|00000600| 6f 66 74 77 61 72 65 2e | 0a 2a 2a 20 0a 2a 2a 20 |oftware.|.** .** |
|00000610| 45 76 65 72 79 6f 6e 65 | 20 69 73 20 67 72 61 6e |Everyone| is gran|
|00000620| 74 65 64 20 70 65 72 6d | 69 73 73 69 6f 6e 20 74 |ted perm|ission t|
|00000630| 6f 20 63 6f 70 79 2c 20 | 6d 6f 64 69 66 79 20 61 |o copy, |modify a|
|00000640| 6e 64 20 72 65 64 69 73 | 74 72 69 62 75 74 65 20 |nd redis|tribute |
|00000650| 74 68 69 73 0a 2a 2a 20 | 4d 65 73 63 68 61 63 68 |this.** |Meschach|
|00000660| 20 4c 69 62 72 61 72 79 | 2c 20 70 72 6f 76 69 64 | Library|, provid|
|00000670| 65 64 3a 0a 2a 2a 20 20 | 31 2e 20 20 41 6c 6c 20 |ed:.** |1. All |
|00000680| 63 6f 70 69 65 73 20 63 | 6f 6e 74 61 69 6e 20 74 |copies c|ontain t|
|00000690| 68 69 73 20 63 6f 70 79 | 72 69 67 68 74 20 6e 6f |his copy|right no|
|000006a0| 74 69 63 65 2e 0a 2a 2a | 20 20 32 2e 20 20 41 6c |tice..**| 2. Al|
|000006b0| 6c 20 6d 6f 64 69 66 69 | 65 64 20 63 6f 70 69 65 |l modifi|ed copie|
|000006c0| 73 20 73 68 61 6c 6c 20 | 63 61 72 72 79 20 61 20 |s shall |carry a |
|000006d0| 6e 6f 74 69 63 65 20 73 | 74 61 74 69 6e 67 20 77 |notice s|tating w|
|000006e0| 68 6f 0a 2a 2a 20 20 20 | 20 20 20 6d 61 64 65 20 |ho.** | made |
|000006f0| 74 68 65 20 6c 61 73 74 | 20 6d 6f 64 69 66 69 63 |the last| modific|
|00000700| 61 74 69 6f 6e 20 61 6e | 64 20 74 68 65 20 64 61 |ation an|d the da|
|00000710| 74 65 20 6f 66 20 73 75 | 63 68 20 6d 6f 64 69 66 |te of su|ch modif|
|00000720| 69 63 61 74 69 6f 6e 2e | 0a 2a 2a 20 20 33 2e 20 |ication.|.** 3. |
|00000730| 20 4e 6f 20 63 68 61 72 | 67 65 20 69 73 20 6d 61 | No char|ge is ma|
|00000740| 64 65 20 66 6f 72 20 74 | 68 69 73 20 73 6f 66 74 |de for t|his soft|
|00000750| 77 61 72 65 20 6f 72 20 | 77 6f 72 6b 73 20 64 65 |ware or |works de|
|00000760| 72 69 76 65 64 20 66 72 | 6f 6d 20 69 74 2e 20 20 |rived fr|om it. |
|00000770| 0a 2a 2a 20 20 20 20 20 | 20 54 68 69 73 20 63 6c |.** | This cl|
|00000780| 61 75 73 65 20 73 68 61 | 6c 6c 20 6e 6f 74 20 62 |ause sha|ll not b|
|00000790| 65 20 63 6f 6e 73 74 72 | 75 65 64 20 61 73 20 63 |e constr|ued as c|
|000007a0| 6f 6e 73 74 72 61 69 6e | 69 6e 67 20 6f 74 68 65 |onstrain|ing othe|
|000007b0| 72 20 73 6f 66 74 77 61 | 72 65 0a 2a 2a 20 20 20 |r softwa|re.** |
|000007c0| 20 20 20 64 69 73 74 72 | 69 62 75 74 65 64 20 6f | distr|ibuted o|
|000007d0| 6e 20 74 68 65 20 73 61 | 6d 65 20 6d 65 64 69 75 |n the sa|me mediu|
|000007e0| 6d 20 61 73 20 74 68 69 | 73 20 73 6f 66 74 77 61 |m as thi|s softwa|
|000007f0| 72 65 2c 20 6e 6f 72 20 | 69 73 20 61 0a 2a 2a 20 |re, nor |is a.** |
|00000800| 20 20 20 20 20 64 69 73 | 74 72 69 62 75 74 69 6f | dis|tributio|
|00000810| 6e 20 66 65 65 20 63 6f | 6e 73 69 64 65 72 65 64 |n fee co|nsidered|
|00000820| 20 61 20 63 68 61 72 67 | 65 2e 0a 2a 2a 0a 2a 2a | a charg|e..**.**|
|00000830| 2a 2a 2a 2a 2a 2a 2a 2a | 2a 2a 2a 2a 2a 2a 2a 2a |********|********|
|00000840| 2a 2a 2a 2a 2a 2a 2a 2a | 2a 2a 2a 2a 2a 2a 2a 2a |********|********|
|00000850| 2a 2a 2a 2a 2a 2a 2a 2a | 2a 2a 2a 2a 2a 2a 2a 2a |********|********|
|00000860| 2a 2a 2a 2a 2a 2a 2a 2a | 2a 2a 2a 2a 2a 2a 2a 2a |********|********|
|00000870| 2a 2a 2a 2a 2a 2a 2a 2a | 2a 2f 0a 0a 0a 2f 2a 20 |********|*/.../* |
|00000880| 69 74 65 72 73 79 6d 2e | 63 20 31 37 2f 30 39 2f |itersym.|c 17/09/|
|00000890| 39 33 20 2a 2f 0a 0a 0a | 2f 2a 20 0a 20 20 49 54 |93 */...|/* . IT|
|000008a0| 45 52 41 54 49 56 45 20 | 4d 45 54 48 4f 44 53 20 |ERATIVE |METHODS |
|000008b0| 2d 20 69 6d 70 6c 65 6d | 65 6e 74 61 74 69 6f 6e |- implem|entation|
|000008c0| 20 6f 66 20 73 65 76 65 | 72 61 6c 20 69 74 65 72 | of seve|ral iter|
|000008d0| 61 74 69 76 65 20 6d 65 | 74 68 6f 64 73 3b 0a 20 |ative me|thods;. |
|000008e0| 20 73 65 65 20 61 6c 73 | 6f 20 69 74 65 72 30 2e | see als|o iter0.|
|000008f0| 63 0a 20 20 2a 2f 0a 0a | 23 69 6e 63 6c 75 64 65 |c. */..|#include|
|00000900| 20 20 20 20 20 20 20 20 | 3c 73 74 64 69 6f 2e 68 | |<stdio.h|
|00000910| 3e 0a 23 69 6e 63 6c 75 | 64 65 09 3c 6d 61 74 68 |>.#inclu|de.<math|
|00000920| 2e 68 3e 0a 23 69 6e 63 | 6c 75 64 65 20 20 20 20 |.h>.#inc|lude |
|00000930| 20 20 20 20 22 6d 61 74 | 72 69 78 2e 68 22 0a 23 | "mat|rix.h".#|
|00000940| 69 6e 63 6c 75 64 65 20 | 20 20 20 20 20 20 20 22 |include | "|
|00000950| 6d 61 74 72 69 78 32 2e | 68 22 0a 23 69 6e 63 6c |matrix2.|h".#incl|
|00000960| 75 64 65 09 22 73 70 61 | 72 73 65 2e 68 22 0a 23 |ude."spa|rse.h".#|
|00000970| 69 6e 63 6c 75 64 65 20 | 20 20 20 20 20 20 20 22 |include | "|
|00000980| 69 74 65 72 2e 68 22 0a | 0a 73 74 61 74 69 63 20 |iter.h".|.static |
|00000990| 63 68 61 72 20 72 63 73 | 69 64 5b 5d 20 3d 20 22 |char rcs|id[] = "|
|000009a0| 24 49 64 3a 20 69 74 65 | 72 73 79 6d 2e 63 2c 76 |$Id: ite|rsym.c,v|
|000009b0| 20 31 2e 31 20 31 39 39 | 34 2f 30 31 2f 31 33 20 | 1.1 199|4/01/13 |
|000009c0| 30 35 3a 33 38 3a 35 39 | 20 64 65 73 20 45 78 70 |05:38:59| des Exp|
|000009d0| 20 24 22 3b 0a 0a 0a 23 | 69 66 64 65 66 20 41 4e | $";...#|ifdef AN|
|000009e0| 53 49 5f 43 0a 56 45 43 | 09 2a 73 70 43 48 73 6f |SI_C.VEC|.*spCHso|
|000009f0| 6c 76 65 28 53 50 4d 41 | 54 20 2a 2c 56 45 43 20 |lve(SPMA|T *,VEC |
|00000a00| 2a 2c 56 45 43 20 2a 29 | 3b 0a 56 45 43 09 2a 74 |*,VEC *)|;.VEC.*t|
|00000a10| 72 69 65 69 67 28 56 45 | 43 20 2a 2c 56 45 43 20 |rieig(VE|C *,VEC |
|00000a20| 2a 2c 4d 41 54 20 2a 29 | 3b 0a 23 65 6c 73 65 0a |*,MAT *)|;.#else.|
|00000a30| 56 45 43 09 2a 73 70 43 | 48 73 6f 6c 76 65 28 29 |VEC.*spC|Hsolve()|
|00000a40| 3b 0a 56 45 43 09 2a 74 | 72 69 65 69 67 28 29 3b |;.VEC.*t|rieig();|
|00000a50| 0a 23 65 6e 64 69 66 0a | 0a 0a 0a 2f 2a 20 69 74 |.#endif.|.../* it|
|00000a60| 65 72 5f 73 70 63 67 20 | 2d 2d 20 61 20 73 69 6d |er_spcg |-- a sim|
|00000a70| 70 6c 65 20 69 6e 74 65 | 72 66 61 63 65 20 74 6f |ple inte|rface to|
|00000a80| 20 69 74 65 72 5f 63 67 | 28 29 20 77 68 69 63 68 | iter_cg|() which|
|00000a90| 20 75 73 65 73 20 73 70 | 61 72 73 65 20 6d 61 74 | uses sp|arse mat|
|00000aa0| 72 69 78 0a 20 20 20 64 | 61 74 61 20 73 74 72 75 |rix. d|ata stru|
|00000ab0| 63 74 75 72 65 73 0a 20 | 20 20 2d 2d 20 61 73 73 |ctures. | -- ass|
|00000ac0| 75 6d 65 73 20 74 68 61 | 74 20 4c 4c 54 20 63 6f |umes tha|t LLT co|
|00000ad0| 6e 74 61 69 6e 73 20 74 | 68 65 20 43 68 6f 6c 65 |ntains t|he Chole|
|00000ae0| 73 6b 79 20 66 61 63 74 | 6f 72 69 73 61 74 69 6f |sky fact|orisatio|
|00000af0| 6e 20 6f 66 20 74 68 65 | 0a 20 20 20 61 63 74 75 |n of the|. actu|
|00000b00| 61 6c 20 70 72 65 63 6f | 6e 64 69 74 69 6f 6e 65 |al preco|nditione|
|00000b10| 72 3b 0a 20 20 20 75 73 | 65 20 61 6c 77 61 79 73 |r;. us|e always|
|00000b20| 20 61 73 20 66 6f 6c 6c | 6f 77 73 3a 0a 20 20 20 | as foll|ows:. |
|00000b30| 78 20 3d 20 69 74 65 72 | 5f 73 70 63 67 28 41 2c |x = iter|_spcg(A,|
|00000b40| 4c 4c 54 2c 62 2c 65 70 | 73 2c 78 2c 6c 69 6d 69 |LLT,b,ep|s,x,limi|
|00000b50| 74 2c 73 74 65 70 73 29 | 3b 0a 20 20 20 6f 72 20 |t,steps)|;. or |
|00000b60| 0a 20 20 20 78 20 3d 20 | 69 74 65 72 5f 73 70 63 |. x = |iter_spc|
|00000b70| 67 28 41 2c 4c 4c 54 2c | 62 2c 65 70 73 2c 56 4e |g(A,LLT,|b,eps,VN|
|00000b80| 55 4c 4c 2c 6c 69 6d 69 | 74 2c 73 74 65 70 73 29 |ULL,limi|t,steps)|
|00000b90| 3b 0a 20 20 20 49 6e 20 | 74 68 65 20 73 65 63 6f |;. In |the seco|
|00000ba0| 6e 64 20 63 61 73 65 20 | 74 68 65 20 73 6f 6c 75 |nd case |the solu|
|00000bb0| 74 69 6f 6e 20 76 65 63 | 74 6f 72 20 69 73 20 63 |tion vec|tor is c|
|00000bc0| 72 65 61 74 65 64 2e 0a | 20 20 20 2a 2f 0a 56 45 |reated..| */.VE|
|00000bd0| 43 20 20 2a 69 74 65 72 | 5f 73 70 63 67 28 41 2c |C *iter|_spcg(A,|
|00000be0| 4c 4c 54 2c 62 2c 65 70 | 73 2c 78 2c 6c 69 6d 69 |LLT,b,ep|s,x,limi|
|00000bf0| 74 2c 73 74 65 70 73 29 | 0a 53 50 4d 41 54 09 2a |t,steps)|.SPMAT.*|
|00000c00| 41 2c 20 2a 4c 4c 54 3b | 0a 56 45 43 09 2a 62 2c |A, *LLT;|.VEC.*b,|
|00000c10| 20 2a 78 3b 0a 64 6f 75 | 62 6c 65 09 65 70 73 3b | *x;.dou|ble.eps;|
|00000c20| 0a 69 6e 74 20 2a 73 74 | 65 70 73 2c 20 6c 69 6d |.int *st|eps, lim|
|00000c30| 69 74 3b 0a 7b 09 0a 20 | 20 20 49 54 45 52 20 2a |it;.{.. | ITER *|
|00000c40| 69 70 3b 0a 20 20 20 0a | 20 20 20 69 70 20 3d 20 |ip;. .| ip = |
|00000c50| 69 74 65 72 5f 67 65 74 | 28 30 2c 30 29 3b 0a 20 |iter_get|(0,0);. |
|00000c60| 20 20 69 70 2d 3e 41 78 | 20 3d 20 28 46 75 6e 5f | ip->Ax| = (Fun_|
|00000c70| 41 78 29 20 73 70 5f 6d | 76 5f 6d 6c 74 3b 0a 20 |Ax) sp_m|v_mlt;. |
|00000c80| 20 20 69 70 2d 3e 41 5f | 70 61 72 20 3d 20 28 76 | ip->A_|par = (v|
|00000c90| 6f 69 64 20 2a 29 41 3b | 0a 20 20 20 69 70 2d 3e |oid *)A;|. ip->|
|00000ca0| 42 78 20 3d 20 28 46 75 | 6e 5f 41 78 29 20 73 70 |Bx = (Fu|n_Ax) sp|
|00000cb0| 43 48 73 6f 6c 76 65 3b | 0a 20 20 20 69 70 2d 3e |CHsolve;|. ip->|
|00000cc0| 42 5f 70 61 72 20 3d 20 | 28 76 6f 69 64 20 2a 29 |B_par = |(void *)|
|00000cd0| 4c 4c 54 3b 0a 20 20 20 | 69 70 2d 3e 69 6e 66 6f |LLT;. |ip->info|
|00000ce0| 20 3d 20 28 46 75 6e 5f | 69 6e 66 6f 29 20 4e 55 | = (Fun_|info) NU|
|00000cf0| 4c 4c 3b 0a 20 20 20 69 | 70 2d 3e 62 20 3d 20 62 |LL;. i|p->b = b|
|00000d00| 3b 0a 20 20 20 69 70 2d | 3e 65 70 73 20 3d 20 65 |;. ip-|>eps = e|
|00000d10| 70 73 3b 0a 20 20 20 69 | 70 2d 3e 6c 69 6d 69 74 |ps;. i|p->limit|
|00000d20| 20 3d 20 6c 69 6d 69 74 | 3b 0a 20 20 20 69 70 2d | = limit|;. ip-|
|00000d30| 3e 78 20 3d 20 78 3b 0a | 20 20 20 69 74 65 72 5f |>x = x;.| iter_|
|00000d40| 63 67 28 69 70 29 3b 0a | 20 20 20 78 20 3d 20 69 |cg(ip);.| x = i|
|00000d50| 70 2d 3e 78 3b 0a 20 20 | 20 69 66 20 28 73 74 65 |p->x;. | if (ste|
|00000d60| 70 73 29 20 2a 73 74 65 | 70 73 20 3d 20 69 70 2d |ps) *ste|ps = ip-|
|00000d70| 3e 73 74 65 70 73 3b 0a | 20 20 20 69 70 2d 3e 73 |>steps;.| ip->s|
|00000d80| 68 61 72 65 64 5f 78 20 | 3d 20 69 70 2d 3e 73 68 |hared_x |= ip->sh|
|00000d90| 61 72 65 64 5f 62 20 3d | 20 54 52 55 45 3b 0a 20 |ared_b =| TRUE;. |
|00000da0| 20 20 69 74 65 72 5f 66 | 72 65 65 28 69 70 29 3b | iter_f|ree(ip);|
|00000db0| 20 20 20 2f 2a 20 72 65 | 6c 65 61 73 65 20 6f 6e | /* re|lease on|
|00000dc0| 6c 79 20 49 54 45 52 20 | 73 74 72 75 63 74 75 72 |ly ITER |structur|
|00000dd0| 65 20 2a 2f 0a 20 20 20 | 72 65 74 75 72 6e 20 78 |e */. |return x|
|00000de0| 3b 09 09 0a 7d 0a 0a 2f | 2a 20 0a 20 20 43 6f 6e |;...}../|* . Con|
|00000df0| 6a 75 67 61 74 65 20 67 | 72 61 64 69 65 6e 74 73 |jugate g|radients|
|00000e00| 20 6d 65 74 68 6f 64 3b | 0a 20 20 2a 2f 0a 56 45 | method;|. */.VE|
|00000e10| 43 20 20 2a 69 74 65 72 | 5f 63 67 28 69 70 29 0a |C *iter|_cg(ip).|
|00000e20| 49 54 45 52 20 2a 69 70 | 3b 0a 7b 0a 20 20 20 73 |ITER *ip|;.{. s|
|00000e30| 74 61 74 69 63 20 56 45 | 43 20 2a 72 20 3d 20 56 |tatic VE|C *r = V|
|00000e40| 4e 55 4c 4c 2c 20 2a 70 | 20 3d 20 56 4e 55 4c 4c |NULL, *p| = VNULL|
|00000e50| 2c 20 2a 71 20 3d 20 56 | 4e 55 4c 4c 2c 20 2a 7a |, *q = V|NULL, *z|
|00000e60| 20 3d 20 56 4e 55 4c 4c | 3b 0a 20 20 20 52 65 61 | = VNULL|;. Rea|
|00000e70| 6c 09 61 6c 70 68 61 2c | 20 62 65 74 61 2c 20 69 |l.alpha,| beta, i|
|00000e80| 6e 6e 65 72 2c 20 6f 6c | 64 5f 69 6e 6e 65 72 2c |nner, ol|d_inner,|
|00000e90| 20 6e 72 65 73 3b 0a 20 | 20 20 56 45 43 20 2a 72 | nres;. | VEC *r|
|00000ea0| 72 3b 20 20 20 2f 2a 20 | 72 72 20 3d 3d 20 72 20 |r; /* |rr == r |
|00000eb0| 6f 72 20 72 72 20 3d 3d | 20 7a 20 2a 2f 0a 20 20 |or rr ==| z */. |
|00000ec0| 20 0a 20 20 20 69 66 20 | 28 69 70 20 3d 3d 20 49 | . if |(ip == I|
|00000ed0| 4e 55 4c 4c 29 0a 20 20 | 20 20 20 65 72 72 6f 72 |NULL). | error|
|00000ee0| 28 45 5f 4e 55 4c 4c 2c | 22 69 74 65 72 5f 63 67 |(E_NULL,|"iter_cg|
|00000ef0| 22 29 3b 0a 20 20 20 69 | 66 20 28 21 69 70 2d 3e |");. i|f (!ip->|
|00000f00| 41 78 20 7c 7c 20 21 69 | 70 2d 3e 62 29 0a 20 20 |Ax || !i|p->b). |
|00000f10| 20 20 20 65 72 72 6f 72 | 28 45 5f 4e 55 4c 4c 2c | error|(E_NULL,|
|00000f20| 22 69 74 65 72 5f 63 67 | 22 29 3b 0a 20 20 20 69 |"iter_cg|");. i|
|00000f30| 66 20 28 20 69 70 2d 3e | 78 20 3d 3d 20 69 70 2d |f ( ip->|x == ip-|
|00000f40| 3e 62 20 29 0a 20 20 20 | 20 20 65 72 72 6f 72 28 |>b ). | error(|
|00000f50| 45 5f 49 4e 53 49 54 55 | 2c 22 69 74 65 72 5f 63 |E_INSITU|,"iter_c|
|00000f60| 67 22 29 3b 0a 20 20 20 | 69 66 20 28 21 69 70 2d |g");. |if (!ip-|
|00000f70| 3e 73 74 6f 70 5f 63 72 | 69 74 29 0a 20 20 20 20 |>stop_cr|it). |
|00000f80| 20 65 72 72 6f 72 28 45 | 5f 4e 55 4c 4c 2c 22 69 | error(E|_NULL,"i|
|00000f90| 74 65 72 5f 63 67 22 29 | 3b 0a 20 20 20 0a 20 20 |ter_cg")|;. . |
|00000fa0| 20 69 66 20 28 20 69 70 | 2d 3e 65 70 73 20 3c 3d | if ( ip|->eps <=|
|00000fb0| 20 30 2e 30 20 29 0a 20 | 20 20 20 20 69 70 2d 3e | 0.0 ). | ip->|
|00000fc0| 65 70 73 20 3d 20 4d 41 | 43 48 45 50 53 3b 0a 20 |eps = MA|CHEPS;. |
|00000fd0| 20 20 0a 20 20 20 72 20 | 3d 20 76 5f 72 65 73 69 | . r |= v_resi|
|00000fe0| 7a 65 28 72 2c 69 70 2d | 3e 62 2d 3e 64 69 6d 29 |ze(r,ip-|>b->dim)|
|00000ff0| 3b 0a 20 20 20 70 20 3d | 20 76 5f 72 65 73 69 7a |;. p =| v_resiz|
|00001000| 65 28 70 2c 69 70 2d 3e | 62 2d 3e 64 69 6d 29 3b |e(p,ip->|b->dim);|
|00001010| 0a 20 20 20 71 20 3d 20 | 76 5f 72 65 73 69 7a 65 |. q = |v_resize|
|00001020| 28 71 2c 69 70 2d 3e 62 | 2d 3e 64 69 6d 29 3b 0a |(q,ip->b|->dim);.|
|00001030| 20 20 20 0a 20 20 20 4d | 45 4d 5f 53 54 41 54 5f | . M|EM_STAT_|
|00001040| 52 45 47 28 72 2c 54 59 | 50 45 5f 56 45 43 29 3b |REG(r,TY|PE_VEC);|
|00001050| 0a 20 20 20 4d 45 4d 5f | 53 54 41 54 5f 52 45 47 |. MEM_|STAT_REG|
|00001060| 28 70 2c 54 59 50 45 5f | 56 45 43 29 3b 0a 20 20 |(p,TYPE_|VEC);. |
|00001070| 20 4d 45 4d 5f 53 54 41 | 54 5f 52 45 47 28 71 2c | MEM_STA|T_REG(q,|
|00001080| 54 59 50 45 5f 56 45 43 | 29 3b 0a 20 20 20 0a 20 |TYPE_VEC|);. . |
|00001090| 20 20 69 66 20 28 69 70 | 2d 3e 42 78 20 21 3d 20 | if (ip|->Bx != |
|000010a0| 28 46 75 6e 5f 41 78 29 | 4e 55 4c 4c 29 20 7b 0a |(Fun_Ax)|NULL) {.|
|000010b0| 20 20 20 20 20 20 7a 20 | 3d 20 76 5f 72 65 73 69 | z |= v_resi|
|000010c0| 7a 65 28 7a 2c 69 70 2d | 3e 62 2d 3e 64 69 6d 29 |ze(z,ip-|>b->dim)|
|000010d0| 3b 0a 20 20 20 20 20 20 | 4d 45 4d 5f 53 54 41 54 |;. |MEM_STAT|
|000010e0| 5f 52 45 47 28 7a 2c 54 | 59 50 45 5f 56 45 43 29 |_REG(z,T|YPE_VEC)|
|000010f0| 3b 0a 20 20 20 20 20 20 | 72 72 20 3d 20 7a 3b 0a |;. |rr = z;.|
|00001100| 20 20 20 7d 0a 20 20 20 | 65 6c 73 65 20 72 72 20 | }. |else rr |
|00001110| 3d 20 72 3b 0a 20 20 20 | 0a 20 20 20 69 66 20 28 |= r;. |. if (|
|00001120| 69 70 2d 3e 78 20 21 3d | 20 56 4e 55 4c 4c 29 20 |ip->x !=| VNULL) |
|00001130| 7b 0a 20 20 20 20 20 20 | 69 66 20 28 69 70 2d 3e |{. |if (ip->|
|00001140| 78 2d 3e 64 69 6d 20 21 | 3d 20 69 70 2d 3e 62 2d |x->dim !|= ip->b-|
|00001150| 3e 64 69 6d 29 0a 09 65 | 72 72 6f 72 28 45 5f 53 |>dim)..e|rror(E_S|
|00001160| 49 5a 45 53 2c 22 69 74 | 65 72 5f 63 67 22 29 3b |IZES,"it|er_cg");|
|00001170| 0a 20 20 20 20 20 20 69 | 70 2d 3e 41 78 28 69 70 |. i|p->Ax(ip|
|00001180| 2d 3e 41 5f 70 61 72 2c | 69 70 2d 3e 78 2c 70 29 |->A_par,|ip->x,p)|
|00001190| 3b 20 20 20 20 09 09 2f | 2a 20 70 20 3d 20 41 2a |; ../|* p = A*|
|000011a0| 78 20 2a 2f 0a 20 20 20 | 20 20 20 76 5f 73 75 62 |x */. | v_sub|
|000011b0| 28 69 70 2d 3e 62 2c 70 | 2c 72 29 3b 09 09 20 09 |(ip->b,p|,r);.. .|
|000011c0| 09 2f 2a 20 72 20 3d 20 | 62 20 2d 20 41 2a 78 20 |./* r = |b - A*x |
|000011d0| 2a 2f 0a 20 20 20 7d 0a | 20 20 20 65 6c 73 65 20 |*/. }.| else |
|000011e0| 7b 20 20 2f 2a 20 69 70 | 2d 3e 78 20 3d 3d 20 30 |{ /* ip|->x == 0|
|000011f0| 20 2a 2f 0a 20 20 20 20 | 20 20 69 70 2d 3e 78 20 | */. | ip->x |
|00001200| 3d 20 76 5f 67 65 74 28 | 69 70 2d 3e 62 2d 3e 64 |= v_get(|ip->b->d|
|00001210| 69 6d 29 3b 0a 20 20 20 | 20 20 20 69 70 2d 3e 73 |im);. | ip->s|
|00001220| 68 61 72 65 64 5f 78 20 | 3d 20 46 41 4c 53 45 3b |hared_x |= FALSE;|
|00001230| 0a 20 20 20 20 20 20 76 | 5f 63 6f 70 79 28 69 70 |. v|_copy(ip|
|00001240| 2d 3e 62 2c 72 29 3b 0a | 20 20 20 7d 0a 20 20 20 |->b,r);.| }. |
|00001250| 0a 20 20 20 6f 6c 64 5f | 69 6e 6e 65 72 20 3d 20 |. old_|inner = |
|00001260| 30 2e 30 3b 0a 20 20 20 | 66 6f 72 20 28 20 69 70 |0.0;. |for ( ip|
|00001270| 2d 3e 73 74 65 70 73 20 | 3d 20 30 3b 20 69 70 2d |->steps |= 0; ip-|
|00001280| 3e 73 74 65 70 73 20 3c | 3d 20 69 70 2d 3e 6c 69 |>steps <|= ip->li|
|00001290| 6d 69 74 3b 20 69 70 2d | 3e 73 74 65 70 73 2b 2b |mit; ip-|>steps++|
|000012a0| 20 29 0a 20 20 20 7b 0a | 20 20 20 20 20 20 69 66 | ). {.| if|
|000012b0| 20 28 20 69 70 2d 3e 42 | 78 20 29 0a 09 28 69 70 | ( ip->B|x )..(ip|
|000012c0| 2d 3e 42 78 29 28 69 70 | 2d 3e 42 5f 70 61 72 2c |->Bx)(ip|->B_par,|
|000012d0| 72 2c 72 72 29 3b 09 09 | 2f 2a 20 72 72 20 3d 20 |r,rr);..|/* rr = |
|000012e0| 42 2a 72 20 2a 2f 0a 20 | 20 20 20 20 20 0a 20 20 |B*r */. | . |
|000012f0| 20 20 20 20 69 6e 6e 65 | 72 20 3d 20 69 6e 5f 70 | inne|r = in_p|
|00001300| 72 6f 64 28 72 72 2c 72 | 29 3b 0a 20 20 20 20 20 |rod(rr,r|);. |
|00001310| 20 6e 72 65 73 20 3d 20 | 73 71 72 74 28 66 61 62 | nres = |sqrt(fab|
|00001320| 73 28 69 6e 6e 65 72 29 | 29 3b 0a 20 20 20 20 20 |s(inner)|);. |
|00001330| 20 69 66 20 28 69 70 2d | 3e 69 6e 66 6f 29 20 69 | if (ip-|>info) i|
|00001340| 70 2d 3e 69 6e 66 6f 28 | 69 70 2c 6e 72 65 73 2c |p->info(|ip,nres,|
|00001350| 72 2c 72 72 29 3b 0a 20 | 20 20 20 20 20 69 66 20 |r,rr);. | if |
|00001360| 28 20 69 70 2d 3e 73 74 | 6f 70 5f 63 72 69 74 28 |( ip->st|op_crit(|
|00001370| 69 70 2c 6e 72 65 73 2c | 72 2c 72 72 29 20 29 20 |ip,nres,|r,rr) ) |
|00001380| 62 72 65 61 6b 3b 0a 20 | 20 20 20 20 20 0a 20 20 |break;. | . |
|00001390| 20 20 20 20 69 66 20 28 | 20 69 70 2d 3e 73 74 65 | if (| ip->ste|
|000013a0| 70 73 20 29 09 2f 2a 20 | 69 66 20 28 20 69 70 2d |ps )./* |if ( ip-|
|000013b0| 3e 73 74 65 70 73 20 3e | 20 30 20 29 20 2e 2e 2e |>steps >| 0 ) ...|
|000013c0| 20 2a 2f 0a 20 20 20 20 | 20 20 7b 0a 09 20 62 65 | */. | {.. be|
|000013d0| 74 61 20 3d 20 69 6e 6e | 65 72 2f 6f 6c 64 5f 69 |ta = inn|er/old_i|
|000013e0| 6e 6e 65 72 3b 0a 09 20 | 70 20 3d 20 76 5f 6d 6c |nner;.. |p = v_ml|
|000013f0| 74 61 64 64 28 72 72 2c | 70 2c 62 65 74 61 2c 70 |tadd(rr,|p,beta,p|
|00001400| 29 3b 0a 20 20 20 20 20 | 20 7d 0a 20 20 20 20 20 |);. | }. |
|00001410| 20 65 6c 73 65 09 09 2f | 2a 20 69 66 20 28 20 69 | else../|* if ( i|
|00001420| 70 2d 3e 73 74 65 70 73 | 20 3d 3d 20 30 20 29 20 |p->steps| == 0 ) |
|00001430| 2e 2e 2e 20 2a 2f 0a 20 | 20 20 20 20 20 7b 0a 09 |... */. | {..|
|00001440| 20 62 65 74 61 20 3d 20 | 30 2e 30 3b 0a 09 20 70 | beta = |0.0;.. p|
|00001450| 20 3d 20 76 5f 63 6f 70 | 79 28 72 72 2c 70 29 3b | = v_cop|y(rr,p);|
|00001460| 0a 09 20 6f 6c 64 5f 69 | 6e 6e 65 72 20 3d 20 30 |.. old_i|nner = 0|
|00001470| 2e 30 3b 0a 20 20 20 20 | 20 20 7d 0a 20 20 20 20 |.0;. | }. |
|00001480| 20 20 28 69 70 2d 3e 41 | 78 29 28 69 70 2d 3e 41 | (ip->A|x)(ip->A|
|00001490| 5f 70 61 72 2c 70 2c 71 | 29 3b 20 20 20 20 20 2f |_par,p,q|); /|
|000014a0| 2a 20 71 20 3d 20 41 2a | 70 20 2a 2f 0a 20 20 20 |* q = A*|p */. |
|000014b0| 20 20 20 61 6c 70 68 61 | 20 3d 20 69 6e 6e 65 72 | alpha| = inner|
|000014c0| 2f 69 6e 5f 70 72 6f 64 | 28 70 2c 71 29 3b 0a 20 |/in_prod|(p,q);. |
|000014d0| 20 20 20 20 20 76 5f 6d | 6c 74 61 64 64 28 69 70 | v_m|ltadd(ip|
|000014e0| 2d 3e 78 2c 70 2c 61 6c | 70 68 61 2c 69 70 2d 3e |->x,p,al|pha,ip->|
|000014f0| 78 29 3b 0a 20 20 20 20 | 20 20 76 5f 6d 6c 74 61 |x);. | v_mlta|
|00001500| 64 64 28 72 2c 71 2c 2d | 61 6c 70 68 61 2c 72 29 |dd(r,q,-|alpha,r)|
|00001510| 3b 0a 20 20 20 20 20 20 | 6f 6c 64 5f 69 6e 6e 65 |;. |old_inne|
|00001520| 72 20 3d 20 69 6e 6e 65 | 72 3b 0a 20 20 20 7d 0a |r = inne|r;. }.|
|00001530| 20 20 20 0a 20 20 20 72 | 65 74 75 72 6e 20 69 70 | . r|eturn ip|
|00001540| 2d 3e 78 3b 0a 7d 0a 0a | 0a 0a 2f 2a 20 69 74 65 |->x;.}..|../* ite|
|00001550| 72 5f 6c 61 6e 63 7a 6f | 73 20 2d 2d 20 72 61 77 |r_lanczo|s -- raw|
|00001560| 20 6c 61 6e 63 7a 6f 73 | 20 61 6c 67 6f 72 69 74 | lanczos| algorit|
|00001570| 68 6d 20 2d 2d 20 6e 6f | 20 72 65 2d 6f 72 74 68 |hm -- no| re-orth|
|00001580| 6f 67 6f 6e 61 6c 69 73 | 61 74 69 6f 6e 0a 20 20 |ogonalis|ation. |
|00001590| 20 2d 2d 20 63 72 65 61 | 74 65 73 20 54 20 6d 61 | -- crea|tes T ma|
|000015a0| 74 72 69 78 20 6f 66 20 | 73 69 7a 65 20 3d 3d 20 |trix of |size == |
|000015b0| 6d 2c 0a 20 20 20 62 75 | 74 20 6e 6f 20 6c 61 72 |m,. bu|t no lar|
|000015c0| 67 65 72 20 74 68 61 6e | 20 62 65 66 6f 72 65 20 |ger than| before |
|000015d0| 62 65 74 61 5f 6b 20 3d | 3d 20 30 0a 20 20 20 2d |beta_k =|= 0. -|
|000015e0| 2d 20 75 73 65 73 20 70 | 61 73 73 65 64 20 72 6f |- uses p|assed ro|
|000015f0| 75 74 69 6e 65 20 74 6f | 20 64 6f 20 6d 61 74 72 |utine to| do matr|
|00001600| 69 78 2d 76 65 63 74 6f | 72 20 6d 75 6c 74 69 70 |ix-vecto|r multip|
|00001610| 6c 69 65 73 20 2a 2f 0a | 76 6f 69 64 09 69 74 65 |lies */.|void.ite|
|00001620| 72 5f 6c 61 6e 63 7a 6f | 73 28 69 70 2c 61 2c 62 |r_lanczo|s(ip,a,b|
|00001630| 2c 62 65 74 61 32 2c 51 | 29 0a 49 54 45 52 20 20 |,beta2,Q|).ITER |
|00001640| 20 20 2a 69 70 3b 0a 56 | 45 43 09 2a 61 2c 20 2a | *ip;.V|EC.*a, *|
|00001650| 62 3b 0a 52 65 61 6c 09 | 2a 62 65 74 61 32 3b 0a |b;.Real.|*beta2;.|
|00001660| 4d 41 54 09 2a 51 3b 0a | 7b 0a 20 20 20 69 6e 74 |MAT.*Q;.|{. int|
|00001670| 09 6a 3b 0a 20 20 20 73 | 74 61 74 69 63 20 56 45 |.j;. s|tatic VE|
|00001680| 43 09 2a 76 20 3d 20 56 | 4e 55 4c 4c 2c 20 2a 77 |C.*v = V|NULL, *w|
|00001690| 20 3d 20 56 4e 55 4c 4c | 2c 20 2a 74 6d 70 20 3d | = VNULL|, *tmp =|
|000016a0| 20 56 4e 55 4c 4c 3b 0a | 20 20 20 52 65 61 6c 09 | VNULL;.| Real.|
|000016b0| 61 6c 70 68 61 2c 20 62 | 65 74 61 2c 20 63 3b 0a |alpha, b|eta, c;.|
|000016c0| 20 20 20 0a 20 20 20 69 | 66 20 28 20 21 20 69 70 | . i|f ( ! ip|
|000016d0| 20 29 0a 20 20 20 20 20 | 65 72 72 6f 72 28 45 5f | ). |error(E_|
|000016e0| 4e 55 4c 4c 2c 22 69 74 | 65 72 5f 6c 61 6e 63 7a |NULL,"it|er_lancz|
|000016f0| 6f 73 22 29 3b 0a 20 20 | 20 69 66 20 28 20 21 20 |os");. | if ( ! |
|00001700| 69 70 2d 3e 41 78 20 7c | 7c 20 21 20 69 70 2d 3e |ip->Ax ||| ! ip->|
|00001710| 78 20 7c 7c 20 21 20 61 | 20 7c 7c 20 21 20 62 20 |x || ! a| || ! b |
|00001720| 29 0a 20 20 20 20 20 65 | 72 72 6f 72 28 45 5f 4e |). e|rror(E_N|
|00001730| 55 4c 4c 2c 22 69 74 65 | 72 5f 6c 61 6e 63 7a 6f |ULL,"ite|r_lanczo|
|00001740| 73 22 29 3b 0a 20 20 20 | 69 66 20 28 20 69 70 2d |s");. |if ( ip-|
|00001750| 3e 6b 20 3c 3d 20 30 20 | 29 0a 20 20 20 20 20 65 |>k <= 0 |). e|
|00001760| 72 72 6f 72 28 45 5f 42 | 4f 55 4e 44 53 2c 22 69 |rror(E_B|OUNDS,"i|
|00001770| 74 65 72 5f 6c 61 6e 63 | 7a 6f 73 22 29 3b 0a 20 |ter_lanc|zos");. |
|00001780| 20 20 69 66 20 28 20 51 | 20 26 26 20 28 20 51 2d | if ( Q| && ( Q-|
|00001790| 3e 6e 20 3c 20 69 70 2d | 3e 78 2d 3e 64 69 6d 20 |>n < ip-|>x->dim |
|000017a0| 7c 7c 20 51 2d 3e 6d 20 | 3c 20 69 70 2d 3e 6b 20 ||| Q->m |< ip->k |
|000017b0| 29 20 29 0a 20 20 20 20 | 20 65 72 72 6f 72 28 45 |) ). | error(E|
|000017c0| 5f 53 49 5a 45 53 2c 22 | 69 74 65 72 5f 6c 61 6e |_SIZES,"|iter_lan|
|000017d0| 63 7a 6f 73 22 29 3b 0a | 20 20 20 0a 20 20 20 61 |czos");.| . a|
|000017e0| 20 3d 20 76 5f 72 65 73 | 69 7a 65 28 61 2c 28 75 | = v_res|ize(a,(u|
|000017f0| 5f 69 6e 74 29 69 70 2d | 3e 6b 29 3b 09 0a 20 20 |_int)ip-|>k);.. |
|00001800| 20 62 20 3d 20 76 5f 72 | 65 73 69 7a 65 28 62 2c | b = v_r|esize(b,|
|00001810| 28 75 5f 69 6e 74 29 28 | 69 70 2d 3e 6b 2d 31 29 |(u_int)(|ip->k-1)|
|00001820| 29 3b 0a 20 20 20 76 20 | 3d 20 76 5f 72 65 73 69 |);. v |= v_resi|
|00001830| 7a 65 28 76 2c 69 70 2d | 3e 78 2d 3e 64 69 6d 29 |ze(v,ip-|>x->dim)|
|00001840| 3b 0a 20 20 20 77 20 3d | 20 76 5f 72 65 73 69 7a |;. w =| v_resiz|
|00001850| 65 28 77 2c 69 70 2d 3e | 78 2d 3e 64 69 6d 29 3b |e(w,ip->|x->dim);|
|00001860| 0a 20 20 20 74 6d 70 20 | 3d 20 76 5f 72 65 73 69 |. tmp |= v_resi|
|00001870| 7a 65 28 74 6d 70 2c 69 | 70 2d 3e 78 2d 3e 64 69 |ze(tmp,i|p->x->di|
|00001880| 6d 29 3b 0a 20 20 20 4d | 45 4d 5f 53 54 41 54 5f |m);. M|EM_STAT_|
|00001890| 52 45 47 28 76 2c 54 59 | 50 45 5f 56 45 43 29 3b |REG(v,TY|PE_VEC);|
|000018a0| 0a 20 20 20 4d 45 4d 5f | 53 54 41 54 5f 52 45 47 |. MEM_|STAT_REG|
|000018b0| 28 77 2c 54 59 50 45 5f | 56 45 43 29 3b 0a 20 20 |(w,TYPE_|VEC);. |
|000018c0| 20 4d 45 4d 5f 53 54 41 | 54 5f 52 45 47 28 74 6d | MEM_STA|T_REG(tm|
|000018d0| 70 2c 54 59 50 45 5f 56 | 45 43 29 3b 0a 20 20 20 |p,TYPE_V|EC);. |
|000018e0| 0a 20 20 20 62 65 74 61 | 20 3d 20 31 2e 30 3b 0a |. beta| = 1.0;.|
|000018f0| 20 20 20 76 5f 7a 65 72 | 6f 28 61 29 3b 0a 20 20 | v_zer|o(a);. |
|00001900| 20 76 5f 7a 65 72 6f 28 | 62 29 3b 0a 20 20 20 69 | v_zero(|b);. i|
|00001910| 66 20 28 51 29 20 6d 5f | 7a 65 72 6f 28 51 29 3b |f (Q) m_|zero(Q);|
|00001920| 0a 20 20 20 0a 20 20 20 | 2f 2a 20 6e 6f 72 6d 61 |. . |/* norma|
|00001930| 6c 69 73 65 20 78 20 61 | 73 20 77 20 2a 2f 0a 20 |lise x a|s w */. |
|00001940| 20 20 63 20 3d 20 76 5f | 6e 6f 72 6d 32 28 69 70 | c = v_|norm2(ip|
|00001950| 2d 3e 78 29 3b 0a 20 20 | 20 69 66 20 28 63 20 3c |->x);. | if (c <|
|00001960| 3d 20 4d 41 43 48 45 50 | 53 29 20 7b 20 2f 2a 20 |= MACHEP|S) { /* |
|00001970| 69 70 2d 3e 78 20 3d 3d | 20 30 20 2a 2f 0a 20 20 |ip->x ==| 0 */. |
|00001980| 20 20 20 20 2a 62 65 74 | 61 32 20 3d 20 30 2e 30 | *bet|a2 = 0.0|
|00001990| 3b 0a 20 20 20 20 20 20 | 72 65 74 75 72 6e 3b 0a |;. |return;.|
|000019a0| 20 20 20 7d 0a 20 20 20 | 65 6c 73 65 20 0a 20 20 | }. |else . |
|000019b0| 20 20 20 73 76 5f 6d 6c | 74 28 31 2e 30 2f 63 2c | sv_ml|t(1.0/c,|
|000019c0| 69 70 2d 3e 78 2c 77 29 | 3b 0a 20 20 20 0a 20 20 |ip->x,w)|;. . |
|000019d0| 20 28 69 70 2d 3e 41 78 | 29 28 69 70 2d 3e 41 5f | (ip->Ax|)(ip->A_|
|000019e0| 70 61 72 2c 77 2c 76 29 | 3b 0a 20 20 20 0a 20 20 |par,w,v)|;. . |
|000019f0| 20 66 6f 72 20 28 20 6a | 20 3d 20 30 3b 20 6a 20 | for ( j| = 0; j |
|00001a00| 3c 20 69 70 2d 3e 6b 3b | 20 6a 2b 2b 20 29 0a 20 |< ip->k;| j++ ). |
|00001a10| 20 20 7b 0a 20 20 20 20 | 20 20 2f 2a 20 73 74 6f | {. | /* sto|
|00001a20| 72 65 20 77 20 69 6e 20 | 51 20 69 66 20 51 20 6e |re w in |Q if Q n|
|00001a30| 6f 74 20 4e 55 4c 4c 20 | 2a 2f 0a 20 20 20 20 20 |ot NULL |*/. |
|00001a40| 20 69 66 20 28 20 51 20 | 29 20 73 65 74 5f 72 6f | if ( Q |) set_ro|
|00001a50| 77 28 51 2c 6a 2c 77 29 | 3b 0a 20 20 20 20 20 20 |w(Q,j,w)|;. |
|00001a60| 0a 20 20 20 20 20 20 61 | 6c 70 68 61 20 3d 20 69 |. a|lpha = i|
|00001a70| 6e 5f 70 72 6f 64 28 77 | 2c 76 29 3b 0a 20 20 20 |n_prod(w|,v);. |
|00001a80| 20 20 20 61 2d 3e 76 65 | 5b 6a 5d 20 3d 20 61 6c | a->ve|[j] = al|
|00001a90| 70 68 61 3b 0a 20 20 20 | 20 20 20 76 5f 6d 6c 74 |pha;. | v_mlt|
|00001aa0| 61 64 64 28 76 2c 77 2c | 2d 61 6c 70 68 61 2c 76 |add(v,w,|-alpha,v|
|00001ab0| 29 3b 0a 20 20 20 20 20 | 20 62 65 74 61 20 3d 20 |);. | beta = |
|00001ac0| 76 5f 6e 6f 72 6d 32 28 | 76 29 3b 0a 20 20 20 20 |v_norm2(|v);. |
|00001ad0| 20 20 69 66 20 28 20 62 | 65 74 61 20 3d 3d 20 30 | if ( b|eta == 0|
|00001ae0| 2e 30 20 29 0a 20 20 20 | 20 20 20 7b 0a 09 20 2a |.0 ). | {.. *|
|00001af0| 62 65 74 61 32 20 3d 20 | 30 2e 30 3b 0a 09 20 72 |beta2 = |0.0;.. r|
|00001b00| 65 74 75 72 6e 3b 0a 20 | 20 20 20 20 20 7d 0a 20 |eturn;. | }. |
|00001b10| 20 20 20 20 20 0a 20 20 | 20 20 20 20 69 66 20 28 | . | if (|
|00001b20| 20 6a 20 3c 20 69 70 2d | 3e 6b 2d 31 20 29 0a 09 | j < ip-|>k-1 )..|
|00001b30| 62 2d 3e 76 65 5b 6a 5d | 20 3d 20 62 65 74 61 3b |b->ve[j]| = beta;|
|00001b40| 0a 20 20 20 20 20 20 76 | 5f 63 6f 70 79 28 77 2c |. v|_copy(w,|
|00001b50| 74 6d 70 29 3b 0a 20 20 | 20 20 20 20 73 76 5f 6d |tmp);. | sv_m|
|00001b60| 6c 74 28 31 2f 62 65 74 | 61 2c 76 2c 77 29 3b 0a |lt(1/bet|a,v,w);.|
|00001b70| 20 20 20 20 20 20 73 76 | 5f 6d 6c 74 28 2d 62 65 | sv|_mlt(-be|
|00001b80| 74 61 2c 74 6d 70 2c 76 | 29 3b 0a 20 20 20 20 20 |ta,tmp,v|);. |
|00001b90| 20 28 69 70 2d 3e 41 78 | 29 28 69 70 2d 3e 41 5f | (ip->Ax|)(ip->A_|
|00001ba0| 70 61 72 2c 77 2c 74 6d | 70 29 3b 0a 20 20 20 20 |par,w,tm|p);. |
|00001bb0| 20 20 76 5f 61 64 64 28 | 76 2c 74 6d 70 2c 76 29 | v_add(|v,tmp,v)|
|00001bc0| 3b 0a 20 20 20 7d 0a 20 | 20 20 2a 62 65 74 61 32 |;. }. | *beta2|
|00001bd0| 20 3d 20 62 65 74 61 3b | 0a 20 20 20 0a 7d 0a 0a | = beta;|. .}..|
|00001be0| 2f 2a 20 69 74 65 72 5f | 73 70 6c 61 6e 63 7a 6f |/* iter_|splanczo|
|00001bf0| 73 20 2d 2d 20 76 65 72 | 73 69 6f 6e 20 74 68 61 |s -- ver|sion tha|
|00001c00| 74 20 75 73 65 73 20 73 | 70 61 72 73 65 20 6d 61 |t uses s|parse ma|
|00001c10| 74 72 69 78 20 64 61 74 | 61 20 73 74 72 75 63 74 |trix dat|a struct|
|00001c20| 75 72 65 20 2a 2f 0a 76 | 6f 69 64 20 20 20 20 69 |ure */.v|oid i|
|00001c30| 74 65 72 5f 73 70 6c 61 | 6e 63 7a 6f 73 28 41 2c |ter_spla|nczos(A,|
|00001c40| 6d 2c 78 30 2c 61 2c 62 | 2c 62 65 74 61 32 2c 51 |m,x0,a,b|,beta2,Q|
|00001c50| 29 0a 53 50 4d 41 54 09 | 2a 41 3b 0a 69 6e 74 20 |).SPMAT.|*A;.int |
|00001c60| 20 20 20 20 6d 3b 0a 56 | 45 43 20 20 20 20 20 2a | m;.V|EC *|
|00001c70| 78 30 2c 20 2a 61 2c 20 | 2a 62 3b 0a 52 65 61 6c |x0, *a, |*b;.Real|
|00001c80| 20 20 20 20 2a 62 65 74 | 61 32 3b 0a 4d 41 54 20 | *bet|a2;.MAT |
|00001c90| 20 20 20 20 2a 51 3b 0a | 7b 09 0a 20 20 20 49 54 | *Q;.|{.. IT|
|00001ca0| 45 52 20 2a 69 70 3b 0a | 20 20 20 0a 20 20 20 69 |ER *ip;.| . i|
|00001cb0| 70 20 3d 20 69 74 65 72 | 5f 67 65 74 28 30 2c 30 |p = iter|_get(0,0|
|00001cc0| 29 3b 0a 20 20 20 69 70 | 2d 3e 73 68 61 72 65 64 |);. ip|->shared|
|00001cd0| 5f 78 20 3d 20 69 70 2d | 3e 73 68 61 72 65 64 5f |_x = ip-|>shared_|
|00001ce0| 62 20 3d 20 54 52 55 45 | 3b 0a 20 20 20 69 70 2d |b = TRUE|;. ip-|
|00001cf0| 3e 41 78 20 3d 20 28 46 | 75 6e 5f 41 78 29 20 73 |>Ax = (F|un_Ax) s|
|00001d00| 70 5f 6d 76 5f 6d 6c 74 | 3b 0a 20 20 20 69 70 2d |p_mv_mlt|;. ip-|
|00001d10| 3e 41 5f 70 61 72 20 3d | 20 28 76 6f 69 64 20 2a |>A_par =| (void *|
|00001d20| 29 20 41 3b 0a 20 20 20 | 69 70 2d 3e 78 20 3d 20 |) A;. |ip->x = |
|00001d30| 78 30 3b 0a 20 20 20 69 | 70 2d 3e 6b 20 3d 20 6d |x0;. i|p->k = m|
|00001d40| 3b 0a 20 20 20 69 74 65 | 72 5f 6c 61 6e 63 7a 6f |;. ite|r_lanczo|
|00001d50| 73 28 69 70 2c 61 2c 62 | 2c 62 65 74 61 32 2c 51 |s(ip,a,b|,beta2,Q|
|00001d60| 29 3b 09 0a 20 20 20 69 | 74 65 72 5f 66 72 65 65 |);.. i|ter_free|
|00001d70| 28 69 70 29 3b 20 20 20 | 2f 2a 20 72 65 6c 65 61 |(ip); |/* relea|
|00001d80| 73 65 20 6f 6e 6c 79 20 | 49 54 45 52 20 73 74 72 |se only |ITER str|
|00001d90| 75 63 74 75 72 65 20 2a | 2f 0a 7d 0a 0a 0a 0a 65 |ucture *|/.}....e|
|00001da0| 78 74 65 72 6e 09 64 6f | 75 62 6c 65 09 66 72 65 |xtern.do|uble.fre|
|00001db0| 78 70 28 29 2c 20 6c 64 | 65 78 70 28 29 3b 0a 0a |xp(), ld|exp();..|
|00001dc0| 2f 2a 20 70 72 6f 64 75 | 63 74 20 2d 2d 20 72 65 |/* produ|ct -- re|
|00001dd0| 74 75 72 6e 73 20 74 68 | 65 20 70 72 6f 64 75 63 |turns th|e produc|
|00001de0| 74 20 6f 66 20 61 20 6c | 6f 6e 67 20 6c 69 73 74 |t of a l|ong list|
|00001df0| 20 6f 66 20 6e 75 6d 62 | 65 72 73 0a 20 20 20 2d | of numb|ers. -|
|00001e00| 2d 20 61 6e 73 77 65 72 | 20 73 74 6f 72 65 64 20 |- answer| stored |
|00001e10| 69 6e 20 6d 61 6e 74 20 | 28 6d 61 6e 74 69 73 73 |in mant |(mantiss|
|00001e20| 61 29 20 61 6e 64 20 65 | 78 70 74 20 28 65 78 70 |a) and e|xpt (exp|
|00001e30| 6f 6e 65 6e 74 | |onent | |
+--------+-------------------------+-------------------------+--------+--------+