function [N,dN]=NURBSDersBasisFunc(ii,pl,u,u_knotl) %nb de derivee dont on a besoin, a priori une seule nders=1; N=zeros(length(u),pl+1); dN=zeros(length(u),pl+1); a=zeros(length(u),nders+1,pl+1); ndu=zeros(length(u),pl+1,pl+1); ders=zeros(length(u),nders+1,pl+1); leftl=zeros(length(u),pl+1); rightl=zeros(length(u),pl+1); N(:,1,1) = 1; ndu(:,1,1)=1; for j = 1:pl leftl(:,j+1) = u - u_knotl(ii+1-j); rightl(:,j+1) = u_knotl(ii+j) - u; saved = 0; for r = 0:(j-1) ndu(:,j+1,r+1) = rightl(:,r+2)+leftl(:,j-r+1); temp = ndu(:,r+1,j)./ndu(:,j+1,r+1); ndu(:,r+1,j+1) = saved + rightl(:,r+2).*temp; saved = leftl(:,j-r+1).*temp; end ndu(:,j+1,j+1) = saved; end for j = 0:pl ders(:,1,j+1)=ndu(:,j+1,pl+1); end for r=0:pl s1=0; s2=1; a(:,1,1)=1; for k=1:nders d=zeros(length(u),1,1); rk=r-k; pk=pl-k; if(r>=k) a(:,s2+1,1)=a(:,s1+1,1)./ndu(:,pk+2,rk+1); d(:,1,1) = a(:,s2+1,1).*ndu(:,rk+1,pk+1); end if(rk>=-1) j1=1; else j1=-rk; end if((r-1)<=pk) j2=k-1; else j2=pl-r; end for j=j1:j2 a(:,s2+1,j+1)=(a(:,s1+1,j+1) - a(:,s1+1,j))./ndu(:,pk+2,rk+j+1); d(:,1,1) = d(:,1,1) + a(:,s2+1,j+1).*ndu(:,rk+j+1,pk+1); end if(r<=pk) a(:,s2+1,k+1) = -a(:,s1+1,k)./ndu(:,pk+2,r+1); d(:,1,1) = d(:,1,1) + a(:,s2+1,k+1).*ndu(:,r+1,pk+1); end ders(:,k+1,r+1)=d(:,1,1); j=s1; s1=s2; s2=j; end end % multiply by correction factor r=pl; for k=1:nders for j=0:pl ders(:,k+1,j+1)=ders(:,k+1,j+1)*r; end r=r*(pl-k); end for j = 0:pl N(:,j+1)=ders(:,1,j+1); dN(:,j+1)=ders(:,2,j+1); end % figure % plot(u,N) % title(sprintf('ii %d',ii+1-pl)) end