Re: R: R: help on Octave

[c.]

On 29 Jun 2012, at 13:08, c. wrote:

> 
> On 29 Jun 2012, at 12:33, Pierpaolo wrote:
> 
>> Hi, 
>> 
>> I created the sparse matrix and sparse array by
>> 
>> s=speye(3*dd);  
>> d=speye(3*dd,1);
>> 
>> This because I don't know how to create an empty sparse matrix:
>> w=s-s;   
>> jj=s-s;
>> fe=d-d;
> 
> to allocate an empty matrix you can do 
> 
> s = sparse (3*dd, 3*dd);
> 
> I see no point in allocating an empty sparse vector.
> 
>> Then I fill the matrix by several for loop like this:
>> 
>> function jac=fja(x)
>> jac=2/(1-x)^3;
>> endfunction
>> 
>> for j=1:(dimytot-4) 
>>      for i=1:(dimxtot-4) 
>>      
>> jj((j+1)*dimxtot+(i+2),(j+1)*dimxtot+(i+2))=feval("fja",vp((j+1)*dimxtot+(i+
>> 2)));
>>      endfor
>> endfor
>> 
>> 
>> [L,U]=lu((w-jj));
>> 
>> err=U\(L\(fe-w*vp));
> 
> It would be _MUCH_ more efficient if you would first create your amtrix in 
> AIJ format and then convert to a sparse matrix, i.e.
> 
> idx = 1;
> for j=1:(dimytot-4) 
>  for i=1:(dimxtot-4) 
>   irow(idx) = (j+1)*dimxtot+(i+2);
>   jcol(idx) = (j+1)*dimxtot+(i+2); 
>   jj(idx++) = feval("fja",vp((j+1)*dimxtot+(i+2)));
>  endfor
> endfor
> 
> jj = sparse (irow, jcol, jj, 3*dd, 3*dd);
> 
> BTW, what's the reason for 
> 
> feval("fja",vp((j+1)*dimxtot+(i+2)))
> 
> why not just 
> 
> fja (vp((j+1)*dimxtot+(i+2)))
> 
> ?
> 
> 
> HTH,
> c.


BTW, it seems you are trying to assemble a large nonlinear system by connecting 
a set of small subsystems, 
like it happens in non-linear Finite Elments problems or in circuit simulation.
The package OCS in Octave-Forge is intended for automating this task, you might 
want to have a look there.

c.