3D グラフィックスのための数学入門 157 ページをコードに．
コードは，kodatuno を参考にさせてもらった．
double	BS_basis	(const double t,c_v_double& knot,const long order,const long i)
{
	if (order <= 0)                  			{	return	0. ;	}
	if (order == 1) {	//	Ni,0(t)
		if (t == knot[knot.size()-1])
		{
			if (knot[i]<=t && t<=knot[i+1])		{	return	1. ;	}
			else                        		{	return	0. ;	}
			}
		else {
			if (knot[i]<=t && t< knot[i+1])		{	return	1. ;	}
			else                        		{	return	0. ;	}
			}
		}
	else
	{
		double	n1 = 0 ;                                	//	(*/0)  -->  n = 0
		double	n2 = 0 ;
		double	d1 = knot[i+order-1] - knot[i  ] ;      	//	T(i+k)  -T(i)
		double	d2 = knot[i+order-0] - knot[i+1] ;      	//	T(i+k+1)-T(i+1)
		if (!::V1_is_near<double<(d1,0))	{
		    	n1 = t - knot[i] ;                      	//	t - T(i)
		    	n1 = n1 * BS_basis(t,knot,order-1,i  ) ;	//	N(i,  k-1)(t)
		    	n1 = n1 / d1 ;
			}
		if (!::V1_is_near<double<(d2,0))	{
		    	n2 = knot[i+order-0] - t ;              	//	T(i+k+1)-t
		    	n2 = n2 * BS_basis(t,knot,order-1,i+1) ;	//	N(i+1,k-1)(t)
		    	n2 = n2 / d2 ;
			}
		return	(n1 + n2) ;
		}
	return	-1. ;
	}
double	BS_basis	(const double t,const long order,const long i)
{
	if (order <= 0)                  			{	return	0. ;	}
	if (order == 1) {	//	Ni,0(t)
		if (i<=t && t< i+1)           			{	return	1. ;	}
		else                            		{	return	0. ;	}
		}
	else
	{
		double	n1 = 0 ;
		double	n2 = 0 ;
		double	d1 = order - 1 ;                    		//	k
		double	d2 = order - 1 ;                    		//	k
		if (!::V1_is_near<double<(d1,0))	{
		    	n1 = t - i ;                        		//	t-i
		    	n1 = n1 * BS_basis(t,order-1,i  ) ; 		//	N(i  ,k-1)
		    	n1 = n1 / d1 ;
			}
		if (!::V1_is_near<double<(d2,0))	{
		    	n2 = (i+order-0) - t ;              		//	i+k+1 - t
		    	n2 = n2 * BS_basis(t,order-1,i+1) ; 		//	N(i+1,k-1)
		    	n2 = n2 / d2 ;
			}
		return	(n1 + n2) ;
		}
	return	-1. ;
	}
	long	pt_count = 7 ;
	v_vd_n.resize(pt_count) ;
	for (size_t index=0 ; index<10 ; index++) {
		for (size_t index_d=0 ;  index_d<div_c ; index_d++) {
			double	t = double(index*div_c + index_d)/div_c ;
			for (size_t i=0 ; i<pt_count ; i++) {
				double	r = ::BS_basis(t,order,i) ;
				Vd3A	v3a = v_vd_n[i] ;
				Vd3	v(v3a.size(),r*100,0) ;
				v3a.push_back(v) ;
				v_vd_n[i] = v3a ;
				}
			}
		}