From mboxrd@z Thu Jan 1 00:00:00 1970 Received: (from majordomo@localhost) by pauillac.inria.fr (8.7.6/8.7.3) id LAA17516; Sat, 1 Feb 2003 11:19:59 +0100 (MET) X-Authentication-Warning: pauillac.inria.fr: majordomo set sender to owner-caml-list@pauillac.inria.fr using -f Received: from nez-perce.inria.fr (nez-perce.inria.fr [192.93.2.78]) by pauillac.inria.fr (8.7.6/8.7.3) with ESMTP id LAA17320 for ; Sat, 1 Feb 2003 11:19:58 +0100 (MET) Received: from shiva.jussieu.fr (shiva.jussieu.fr [134.157.0.129]) by nez-perce.inria.fr (8.11.1/8.11.1) with ESMTP id h11AJvP11351 for ; Sat, 1 Feb 2003 11:19:57 +0100 (MET) Received: from ibm3.cicrp.jussieu.fr (ibm3.cicrp.jussieu.fr [134.157.15.3]) by shiva.jussieu.fr (8.12.5/jtpda-5.4) with ESMTP id h11AJvYY098146 ; Sat, 1 Feb 2003 11:19:57 +0100 (CET) Received: from ibm1.cicrp.jussieu.fr (ibm1.cicrp.jussieu.fr [134.157.15.1]) by ibm3.cicrp.jussieu.fr (8.8.8/jtpda/mob-V8) with ESMTP id LAA21304 ; Sat, 1 Feb 2003 11:18:49 +0100 Received: from localhost (fernande@localhost) by ibm1.cicrp.jussieu.fr (8.8.8/jtpda/mob-v8) with ESMTP id LAA1048694 ; Sat, 1 Feb 2003 11:18:48 +0100 Date: Sat, 1 Feb 2003 11:18:48 +0100 (NFT) From: Diego Olivier Fernandez Pons To: Brian Hurt cc: caml-list@inria.fr Subject: Linear systems (was Re: [Caml-list] @, List.append, and tail recursion) In-Reply-To: Message-ID: MIME-Version: 1.0 Content-Type: TEXT/PLAIN; charset=iso-8859-1 Content-Transfer-Encoding: QUOTED-PRINTABLE X-Antivirus: scanned by sophie at shiva.jussieu.fr Sender: owner-caml-list@pauillac.inria.fr Precedence: bulk Bonjour, Sorry for the buggy code, anyway you should be able to fix it. > In production I wouldn't be surprised to see systems of 30,000+ > equations in 30,000+ variables. The Jacobian matrix is going to be > very sparse- I expect the average row to have 20 or fewer non-zero > elements. Thus the attraction to sparse vectors. I'm going to be > solving it via Gaussian elimination (the Jacobian is likely to be > malformed in multiple ways, meaning I can't use any iterative method > I know of. And yes, I've looked at iterative methods as advanced as > GMRES- they don't work). I think that in general I can bound the > size of vectors I'm producing. But there are degenerate cases where > I could get above 30K non-zero elements in a vector. A 30 000 x 30 000 system of equations is not a toy program. Then, do not expect to solve it straightforwardly. You will obviously hit all system limits (stack, cache, ...). You have to understand that if you really want a robuts program in all cases, you will have to work : You will have to design several data structures for every operation you want (to keep both efficiency and generality) and transformation functions. That is the case for example in the Gr=F6bner basis system I pointed out : several monomial orderings with many transformation functions ... You will have to use all Caml properties (functional, imperative ...) You also need to understand roughtly how does the Caml compiler work (boxing/ unboxing, garbage collection, ...) > But I want the rare case to *work* correctly, even if inefficiently. Wit= h > the "naive" non-tail-recursive implementation doesn't. Somewhere above > 32K elements in the list the recursion trips the stack overflow. Change > the problem in some minor way, and suddenly I'm not generating lists with > 32K elements in them, maybe just 30K elements in them, and everything > works OK. Your main problem is that sometimes your lists may be too big for the data structure you have chosen. Then, the best thing to do is to have two data structures, one for small (most of all) systems, a second one for huge (but rare) systems. - lists for small sparse vectors - (say) hashtables for huge sparse vectors You just need to write you own hashtable data structure (because you can then profit of Caml's float array unboxing optimisation by separate collision resolution) type size =3D int type vector =3D | List of size * (int * float) list | Hashtable of myhashtbl Most of your code will be purely functional and the program will not break on large data. One more point : floating arithmetic may produce incorrect value by error propagation, even for small systems. If you do want you system to work properly for all data (even not well scaled), you will have to use specific algorithms (pivot choice rules, error bounding, ...). Then, list representation may not be apropriate. > O(1), or O(log n)? Most tree operations are O(log n). It is easy to design a data structure with O(log n) acces to any element and O(1) acces to the first one : imagine a list of increasing perfect trees of size 2^k Diego Olivier ------------------- To unsubscribe, mail caml-list-request@inria.fr Archives: http://caml.inria.fr Bug reports: http://caml.inria.fr/bin/caml-bugs FAQ: http://caml.inria.fr/FAQ/ Beginner's list: http://groups.yahoo.com/group/ocaml_beginners