From mboxrd@z Thu Jan 1 00:00:00 1970 Return-Path: Received: (qmail 18631 invoked by alias); 6 Oct 2009 19:59:49 -0000 Received: (qmail 18623 invoked by uid 22791); 6 Oct 2009 19:59:48 -0000 X-SWARE-Spam-Status: No, hits=-10.1 required=5.0 tests=AWL,BAYES_00,RCVD_IN_DNSWL_HI,SPF_PASS X-Spam-Check-By: sourceware.org Received: from proofpoint1.lanl.gov (HELO proofpoint1.lanl.gov) (204.121.3.25) by sourceware.org (qpsmtpd/0.43rc1) with ESMTP; Tue, 06 Oct 2009 19:59:43 +0000 Received: from mailrelay1.lanl.gov (mailrelay1.lanl.gov [128.165.4.101]) by proofpoint1.lanl.gov (8.14.3/8.14.3) with ESMTP id n96JxeTW001425 for ; Tue, 6 Oct 2009 13:59:41 -0600 Received: from alvie-mail.lanl.gov (alvie-mail.lanl.gov [128.165.4.110]) by mailrelay1.lanl.gov (Postfix) with ESMTP id E597418010C for ; Tue, 6 Oct 2009 13:59:40 -0600 (MDT) Received: from localhost (localhost.localdomain [127.0.0.1]) by alvie-mail.lanl.gov (Postfix) with ESMTP id E474A7D00D9 for ; Tue, 6 Oct 2009 13:59:40 -0600 (MDT) X-NIE-2-Virus-Scanner: amavisd-new at alvie-mail.lanl.gov Received: from [130.55.124.157] (manticore.lanl.gov [130.55.124.157]) by alvie-mail.lanl.gov (Postfix) with ESMTP id CBBEB7D00D8 for ; Tue, 6 Oct 2009 13:59:40 -0600 (MDT) Subject: Re: containers tentative design summary From: Gerard Jungman To: gsl-discuss@sourceware.org In-Reply-To: <7f1eaee30910051645v5654353fj53416ed88b98844@mail.gmail.com> References: <1254708349.18519.4.camel@ForbiddenPlanet> <7f1eaee30910050750l738876b1p41e6bd8ae5aa6d16@mail.gmail.com> <1254783367.28192.98.camel@manticore.lanl.gov> <7f1eaee30910051645v5654353fj53416ed88b98844@mail.gmail.com> Content-Type: text/plain Date: Tue, 06 Oct 2009 19:59:00 -0000 Message-Id: <1254859325.28192.209.camel@manticore.lanl.gov> Mime-Version: 1.0 Content-Transfer-Encoding: 7bit X-Proofpoint-Virus-Version: vendor=fsecure engine=1.12.8161:2.4.5,1.2.40,4.0.166 definitions=2009-10-06_17:2009-09-29,2009-10-06,2009-10-06 signatures=0 Mailing-List: contact gsl-discuss-help@sourceware.org; run by ezmlm Precedence: bulk List-Id: List-Subscribe: List-Archive: List-Post: List-Help: , Sender: gsl-discuss-owner@sourceware.org X-SW-Source: 2009-q4/txt/msg00007.txt.bz2 On Mon, 2009-10-05 at 19:45 -0400, James Bergstra wrote: > I'm thinking of something like > > struct marray { int rank; int * dims; int * strides; dtype * base; }; > marray_add_inplace( const marray a, const marray b, marray c); //N.B. > the passing on the stack > > Inside that function, we can check for special cases of interest and > optimize them as necessary. This is not so clear. You can't have a test inside the indexing functions; that would obviously kill you. So we still need a way to pick the correct indexing function at compile-time; I think this means that we need the fixed-rank types. And the indexing calculation is best done with 'dims' and 'strides' on the stack, rather than through the pointer indirection. Maybe you are suggesting the above as an implementation detail, to be wrapped by the fixed-rank types? That would allow for compile-time selection of the indexing, but it would still probably force an indirection through those pointers. > One disadvantage of such a structure I think you'll agree, is the > internal pointers, which interfere with convenient stack allocation. Yes. I rejected internal pointers for that reason. > A somewhat ugly alternative (which I am nonetheless suggesting) would > be the following... > > { int rank; dtype * base; int dim0, dim1, dim2; int str0, str1, str2; > int * extra_dims; int * extra_strides; }; > > For low-rank multiarrays, this structure can be manipulated on the > stack, so it is possible to write expressions like > > marray_add_inplace(marray_reshape2(a, 5, 50), marray_slice(b, 17), > marray_transpose(c)); // add a reshaped view of a to the 17th > matrix-slice of b, storing result in a transposed view of c. > > At the same time, arbitrarily high ranks are supported if the user is > willing to put up with the syntax of heap-based manipulation. This is closer to what I have done, although I just separated the fixed-rank and the arbitrary rank (heap-based) cases. In the end, that seemed like the cleanest choice. > (Of course we can haggle over how many extra ranks can be built into > the static part of the structure.) The way I have done it, I avoid this question. Here are some snippets from what I have, just so you get an idea of how I am thinking about it. I'll work downwards, from the multi-arrays, down to the implementation details. First, we must have separate types for the fixed rank cases, in order to insure that we select the correct indexing function at compile-time. So... typedef struct { ... } gsl_marray_1; // rank 1 typedef struct { ... } gsl_marray_2; // rank 2 typedef struct { ... } gal_marray_3; // rank 3 etc. The basic ingredient for indexing is a dimension and stride in an index place. So... typedef struct { size_t dimens; // dimension in this index place size_t stride; // addressing stride for this index } gsl_marray_idx; Given this, the marray declarations can look like: typedef struct { int rank; double * data; gsl_marray_idx idx[1]; ... } gsl_marray_1 typedef struct { int rank; double * data; gsl_marray_idx idx[2]; ... } gsl_marray_2 etc. Now, at the level of implementations, we write everything in terms of (rank, gsl_marray_idx[]), avoiding an intermediate type. We can freely mix generic versions of functions and specializations. Generic versions handle the general case, for operations which are not performance critical (like formatted i/o and stuff...). Specializations handle the things that matter, like indexing. Examples: // some "low-level" io implementation function; // handles the access generically gsl_marray_fread_raw(FILE * f, double * data, int rank, const gsl_marray_idx idx[]); // indexing (offset calculation) requires specialization inline size_t gsl_marray_idx_1_offset(const gsl_marray_idx idx[], const int index[]) { return index[0]*idx[0].stride; } inline size_t gsl_marray_idx_2_offset(const gsl_marray_idx idx[], const int index[]) { return index[0]*idx[0].stride + index[1]*idx[1].stride; } etc. inline gsl_marray_1_offset(const gsl_marray_1 * ma, const int index[]) { return gsl_marray_idx_1_offset(ma->idx, index); } You get the idea. These are taken (almost) verbatim from the code that I have written so far. Finally, for the generic case we have something like: typedef struct { int rank; gsl_marray_idx * idx; double * data; ... }; with the obvious heap allocating constructor, etc. The same low-level generic functions as above are also used to implement the class-specific functions here. Anyway, the whole thing seems to be layered properly, and there are no "tricky" bits. No rocket science. Just following my nose. I will shortly produce a complete package for people to examine. But you know how it is. As David Byrne once sang, "you can't see it 'til it's finished...". In the meantime, keep suggesting stuff. Discussion always helps. -- G. Jungman