Skip to main content
replaced http://mathematica.stackexchange.com/ with https://mathematica.stackexchange.com/
Source Link

We're benchmarking both time and memory usage using my timeAndMemoryBenchmarkPlots functiontimeAndMemoryBenchmarkPlots function.

Since lazy lists cache their values and benchmark measures running time multiple times, we clear cache after each usage, using LazyListBlock (as per comment by Leonidcomment by Leonid).

We're benchmarking both time and memory usage using my timeAndMemoryBenchmarkPlots function.

Since lazy lists cache their values and benchmark measures running time multiple times, we clear cache after each usage, using LazyListBlock (as per comment by Leonid).

We're benchmarking both time and memory usage using my timeAndMemoryBenchmarkPlots function.

Since lazy lists cache their values and benchmark measures running time multiple times, we clear cache after each usage, using LazyListBlock (as per comment by Leonid).

Bounty Awarded with 100 reputation awarded by ciao
Replace `Module` with `With` where possible. Use `LazyListBlock`. Add benchmarks for lazy lists with `k=8`.
Source Link
jkuczm
  • 15.2k
  • 2
  • 55
  • 89
Needs["Streaming`"] ClearAll[lazyWeakCompositions] lazyWeakCompositions[ n_Integer?NonNegative, 0, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[] ] := LazyListCreate[{}, chunkSize, opts] lazyWeakCompositions[ n_Integer?NonNegative, 1, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[] ] := LazyListCreate[{{n}}, chunkSize, opts] lazyWeakCompositions[ n_Integer?NonNegative, k_Integer /; k >= 2, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[] ] := With[{length = Binomial[n + k - 1, n]}, Module[ { active = False, compositionsLeft = length, prevComposition = Join[ConstantArray[0, k - 2], {-1, n + 1}] // Developer`ToPackedArray, lastNonZeroIndex = {k - 1} // Developer`ToPackedArray } , LazyListCreate[ IteratorCreate[ ListIterator, (active = True) &, Module[With[{realChunkSize = Min[chunkSize, compositionsLeft]}, If[realChunkSize === 0, {} ,(* else *) compositionsLeft -= realChunkSize; nextWeakCompositionsChunk[ n, k, prevComposition, lastNonZeroIndex, realChunkSize ] ] ] &, TrueQ[active] &, Remove[ active, compositionsLeft, prevComposition, lastNonZeroIndex ] & ], chunkSize, opts, "Length" -> length, "FiniteList" -> True ] ] ] 
  1. All compositions generated at once:

     ClearAll[scanWeakCompositionsPermPartNonJoined] scanWeakCompositionsPermPartNonJoined[n_, k_] := ( Scan[IdentityLazyListBlock@Scan[Identity, weakCompositionsPermPartNonJoined[n, k], {2}]; DeleteFile@FileNames["*.mx", $StreamingCacheBase] )] 
  2. Compositions generated by permuting lazy partitions. Partitions are taken one by one, so in single chunk we have all permutations of given partition.

     ClearAll[lazyWeakCompositionsPermPart, scanLazyWeakCompositionsPermPart] lazyWeakCompositionsPermPart[n_, k_] := Permutations /@ lazyIntegerPartitions[n, {k}, Range[0, n], 1] scanLazyWeakCompositionsPermPart[n_, k_] := ( Scan[Scan[Identity]LazyListBlock@Scan[Scan[Identity], lazyWeakCompositionsPermPart[n, k]];  DeleteFile@FileNames["*.mx", $StreamingCacheBase] )k]] 
  3. Compositions generated using our lazyWeakCompositions function in chunks of default length: 10^5.

     ClearAll[scanLazyWeakCompositions] scanLazyWeakCompositions[n_, k_]k_, chunkSize_:= (10^5] := Scan[IdentityLazyListBlock@Scan[Identity, lazyWeakCompositions[n, k]];  DeleteFile@FileNames["*.mx"k, $StreamingCacheBase] )chunkSize]] 

Since lazy lists cache their values and benchmark measures running time multiple times, we clear cache after each usage, using LazyListBlock (as per comment by Leonid).

$functionsList = {scanWeakCompositionsPermPartNonJoined, scanLazyWeakCompositionsPermPart, scanLazyWeakCompositions}; timeAndMemoryBenchmarkPlots[fixedKFunctionsList[8], Identity, {1, 2, 3, 4, 5, 6, 7, 8, 9, 11, 12, 14, 16, 18, 21, 24, 28, 32, 37, 42}, TimeConstraint -> Infinity, MemoryConstraint -> 2^31] timeAndMemoryBenchmarkPlots[fixedKFunctionsList[16], Identity, Range[16], TimeConstraint -> Infinity, MemoryConstraint -> 2^31] timeAndMemoryBenchmarkPlots[fixedKFunctionsList[32], Identity, Range[8], TimeConstraint -> Infinity, MemoryConstraint -> 2^31] 

enter image description herebenchmarks lazyWeakCompositions fixed k

For k=8 number of permutations per partition is relatively small, so scanLazyWeakCompositionsPermPart uses a lot of relatively small chunks, in tested range of n. That's why it's considerably slower and it's memory usage increases slowly with n.

ClearAll[nextWeakCompositionsChunk] nextWeakCompositionsChunk = Compile[{{n, _Integer}, {k, _Integer}, {prevComposition, _Integer, 1}, {chunkSize, _Integer}}, Module[{lastNonZeroPosition = k, composition = prevComposition}, While[composition[[lastNonZeroPosition]] == 0, --lastNonZeroPosition]; Table[ ++composition[[lastNonZeroPosition - 1]]; If[lastNonZeroPosition == k, --composition[[-1]] ,(* else *) composition[[-1]] = composition[[lastNonZeroPosition]] - 1; composition[[lastNonZeroPosition]] = 0 ]; If[composition[[-1]] == 0, --lastNonZeroPosition ,(* else *) lastNonZeroPosition = k ]; composition , chunkSize ] ], RuntimeOptions -> "Speed", CompilationTarget -> "WVM" ]; ClearAll[weakCompositions] weakCompositions[n_Integer?NonNegative, 0] := {} weakCompositions[n_Integer?NonNegative, 1] := {{n}} weakCompositions[n_Integer?NonNegative, k_Integer /; k >= 2] := nextWeakCompositionsChunk[n, k, Join[ConstantArray[0, k - 2], {-1, n + 1}], Binomial[n + k - 1, n]] Needs["Streaming`"] ClearAll[lazyWeakCompositions] lazyWeakCompositions[n_Integer?NonNegative, 0, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[]] := LazyListCreate[{}, chunkSize, opts] lazyWeakCompositions[n_Integer?NonNegative, 1, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[]] := LazyListCreate[{{n}}, chunkSize, opts] lazyWeakCompositions[n_Integer?NonNegative, k_Integer /; k >= 2, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[]] := With[{length = Binomial[n + k - 1, n]}, Module[{active = False, compositionsLeft = length, prevComposition = Join[ConstantArray[0, k - 2], {-1, n + 1}]}, LazyListCreate[ IteratorCreate[ ListIterator, (active = True) &, Module[With[{realChunkSize = Min[chunkSize, compositionsLeft], taken}, If[realChunkSize === 0, {} ,(* else **else*) compositionsLeft -= realChunkSize; With[{taken = nextWeakCompositionsChunk[n, k, prevComposition, realChunkSize];realChunkSize]},   prevComposition = Last[taken];   taken ] ] ] &, TrueQ[active] &, Remove[active, compositionsLeft, prevComposition] & ], chunkSize, opts, "Length" -> length, "FiniteList" -> True ] ] ] 
Needs["Streaming`"] ClearAll[lazyWeakCompositions] lazyWeakCompositions[ n_Integer?NonNegative, 0, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[] ] := LazyListCreate[{}, chunkSize, opts] lazyWeakCompositions[ n_Integer?NonNegative, 1, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[] ] := LazyListCreate[{{n}}, chunkSize, opts] lazyWeakCompositions[ n_Integer?NonNegative, k_Integer /; k >= 2, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[] ] := With[{length = Binomial[n + k - 1, n]}, Module[ { active = False, compositionsLeft = length, prevComposition = Join[ConstantArray[0, k - 2], {-1, n + 1}] // Developer`ToPackedArray, lastNonZeroIndex = {k - 1} // Developer`ToPackedArray } , LazyListCreate[ IteratorCreate[ ListIterator, (active = True) &, Module[{realChunkSize = Min[chunkSize, compositionsLeft]}, If[realChunkSize === 0, {} ,(* else *) compositionsLeft -= realChunkSize; nextWeakCompositionsChunk[ n, k, prevComposition, lastNonZeroIndex, realChunkSize ] ] ] &, TrueQ[active] &, Remove[ active, compositionsLeft, prevComposition, lastNonZeroIndex ] & ], chunkSize, opts, "Length" -> length, "FiniteList" -> True ] ] ] 
  1. All compositions generated at once:

     ClearAll[scanWeakCompositionsPermPartNonJoined] scanWeakCompositionsPermPartNonJoined[n_, k_] := ( Scan[Identity, weakCompositionsPermPartNonJoined[n, k], {2}]; DeleteFile@FileNames["*.mx", $StreamingCacheBase] ) 
  2. Compositions generated by permuting lazy partitions. Partitions are taken one by one, so in single chunk we have all permutations of given partition.

     ClearAll[lazyWeakCompositionsPermPart, scanLazyWeakCompositionsPermPart] lazyWeakCompositionsPermPart[n_, k_] := Permutations /@ lazyIntegerPartitions[n, {k}, Range[0, n], 1] scanLazyWeakCompositionsPermPart[n_, k_] := ( Scan[Scan[Identity], lazyWeakCompositionsPermPart[n, k]];  DeleteFile@FileNames["*.mx", $StreamingCacheBase] ) 
  3. Compositions generated using our lazyWeakCompositions function in chunks of default length: 10^5.

     ClearAll[scanLazyWeakCompositions] scanLazyWeakCompositions[n_, k_] := ( Scan[Identity, lazyWeakCompositions[n, k]];  DeleteFile@FileNames["*.mx", $StreamingCacheBase] ) 

Since lazy lists cache their values and benchmark measures running time multiple times, we clear cache after each usage.

timeAndMemoryBenchmarkPlots[fixedKFunctionsList[16], Identity, Range[16], TimeConstraint -> Infinity, MemoryConstraint -> 2^31] timeAndMemoryBenchmarkPlots[fixedKFunctionsList[32], Identity, Range[8], TimeConstraint -> Infinity, MemoryConstraint -> 2^31] 

enter image description here

ClearAll[nextWeakCompositionsChunk] nextWeakCompositionsChunk = Compile[{{n, _Integer}, {k, _Integer}, {prevComposition, _Integer, 1}, {chunkSize, _Integer}}, Module[{lastNonZeroPosition = k, composition = prevComposition}, While[composition[[lastNonZeroPosition]] == 0, --lastNonZeroPosition]; Table[ ++composition[[lastNonZeroPosition - 1]]; If[lastNonZeroPosition == k, --composition[[-1]] ,(* else *) composition[[-1]] = composition[[lastNonZeroPosition]] - 1; composition[[lastNonZeroPosition]] = 0 ]; If[composition[[-1]] == 0, --lastNonZeroPosition ,(* else *) lastNonZeroPosition = k ]; composition , chunkSize ] ], RuntimeOptions -> "Speed", CompilationTarget -> "WVM" ]; ClearAll[weakCompositions] weakCompositions[n_Integer?NonNegative, 0] := {} weakCompositions[n_Integer?NonNegative, 1] := {{n}} weakCompositions[n_Integer?NonNegative, k_Integer /; k >= 2] := nextWeakCompositionsChunk[n, k, Join[ConstantArray[0, k - 2], {-1, n + 1}], Binomial[n + k - 1, n]] Needs["Streaming`"] ClearAll[lazyWeakCompositions] lazyWeakCompositions[n_Integer?NonNegative, 0, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[]] := LazyListCreate[{}, chunkSize, opts] lazyWeakCompositions[n_Integer?NonNegative, 1, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[]] := LazyListCreate[{{n}}, chunkSize, opts] lazyWeakCompositions[n_Integer?NonNegative, k_Integer /; k >= 2, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[]] := With[{length = Binomial[n + k - 1, n]}, Module[{active = False, compositionsLeft = length, prevComposition = Join[ConstantArray[0, k - 2], {-1, n + 1}]}, LazyListCreate[ IteratorCreate[ ListIterator, (active = True) &, Module[{realChunkSize = Min[chunkSize, compositionsLeft], taken}, If[realChunkSize === 0, {} ,(* else *) compositionsLeft -= realChunkSize; taken = nextWeakCompositionsChunk[n, k, prevComposition, realChunkSize]; prevComposition = Last[taken]; taken ] ] &, TrueQ[active] &, Remove[active, compositionsLeft, prevComposition] & ], chunkSize, opts, "Length" -> length, "FiniteList" -> True ] ] ] 
Needs["Streaming`"] ClearAll[lazyWeakCompositions] lazyWeakCompositions[ n_Integer?NonNegative, 0, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[] ] := LazyListCreate[{}, chunkSize, opts] lazyWeakCompositions[ n_Integer?NonNegative, 1, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[] ] := LazyListCreate[{{n}}, chunkSize, opts] lazyWeakCompositions[ n_Integer?NonNegative, k_Integer /; k >= 2, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[] ] := With[{length = Binomial[n + k - 1, n]}, Module[ { active = False, compositionsLeft = length, prevComposition = Join[ConstantArray[0, k - 2], {-1, n + 1}] // Developer`ToPackedArray, lastNonZeroIndex = {k - 1} // Developer`ToPackedArray } , LazyListCreate[ IteratorCreate[ ListIterator, (active = True) &, With[{realChunkSize = Min[chunkSize, compositionsLeft]}, If[realChunkSize === 0, {} ,(* else *) compositionsLeft -= realChunkSize; nextWeakCompositionsChunk[ n, k, prevComposition, lastNonZeroIndex, realChunkSize ] ] ] &, TrueQ[active] &, Remove[ active, compositionsLeft, prevComposition, lastNonZeroIndex ] & ], chunkSize, opts, "Length" -> length, "FiniteList" -> True ] ] ] 
  1. All compositions generated at once:

     ClearAll[scanWeakCompositionsPermPartNonJoined] scanWeakCompositionsPermPartNonJoined[n_, k_] := LazyListBlock@Scan[Identity, weakCompositionsPermPartNonJoined[n, k], {2}] 
  2. Compositions generated by permuting lazy partitions. Partitions are taken one by one, so in single chunk we have all permutations of given partition.

     ClearAll[lazyWeakCompositionsPermPart, scanLazyWeakCompositionsPermPart] lazyWeakCompositionsPermPart[n_, k_] := Permutations /@ lazyIntegerPartitions[n, {k}, Range[0, n], 1] scanLazyWeakCompositionsPermPart[n_, k_] := LazyListBlock@Scan[Scan[Identity], lazyWeakCompositionsPermPart[n, k]] 
  3. Compositions generated using our lazyWeakCompositions function in chunks of default length: 10^5.

     ClearAll[scanLazyWeakCompositions] scanLazyWeakCompositions[n_, k_, chunkSize_: 10^5] := LazyListBlock@Scan[Identity, lazyWeakCompositions[n, k, chunkSize]] 

Since lazy lists cache their values and benchmark measures running time multiple times, we clear cache after each usage, using LazyListBlock (as per comment by Leonid).

$functionsList = {scanWeakCompositionsPermPartNonJoined, scanLazyWeakCompositionsPermPart, scanLazyWeakCompositions}; timeAndMemoryBenchmarkPlots[fixedKFunctionsList[8], Identity, {1, 2, 3, 4, 5, 6, 7, 8, 9, 11, 12, 14, 16, 18, 21, 24, 28, 32, 37, 42}, TimeConstraint -> Infinity, MemoryConstraint -> 2^31] timeAndMemoryBenchmarkPlots[fixedKFunctionsList[16], Identity, Range[16], TimeConstraint -> Infinity, MemoryConstraint -> 2^31] timeAndMemoryBenchmarkPlots[fixedKFunctionsList[32], Identity, Range[8], TimeConstraint -> Infinity, MemoryConstraint -> 2^31] 

benchmarks lazyWeakCompositions fixed k

For k=8 number of permutations per partition is relatively small, so scanLazyWeakCompositionsPermPart uses a lot of relatively small chunks, in tested range of n. That's why it's considerably slower and it's memory usage increases slowly with n.

ClearAll[nextWeakCompositionsChunk] nextWeakCompositionsChunk = Compile[{{n, _Integer}, {k, _Integer}, {prevComposition, _Integer, 1}, {chunkSize, _Integer}}, Module[{lastNonZeroPosition = k, composition = prevComposition}, While[composition[[lastNonZeroPosition]] == 0, --lastNonZeroPosition]; Table[ ++composition[[lastNonZeroPosition - 1]]; If[lastNonZeroPosition == k, --composition[[-1]] ,(* else *) composition[[-1]] = composition[[lastNonZeroPosition]] - 1; composition[[lastNonZeroPosition]] = 0 ]; If[composition[[-1]] == 0, --lastNonZeroPosition ,(* else *) lastNonZeroPosition = k ]; composition , chunkSize ] ], RuntimeOptions -> "Speed", CompilationTarget -> "WVM" ]; ClearAll[weakCompositions] weakCompositions[n_Integer?NonNegative, 0] := {} weakCompositions[n_Integer?NonNegative, 1] := {{n}} weakCompositions[n_Integer?NonNegative, k_Integer /; k >= 2] := nextWeakCompositionsChunk[n, k, Join[ConstantArray[0, k - 2], {-1, n + 1}], Binomial[n + k - 1, n]] Needs["Streaming`"] ClearAll[lazyWeakCompositions] lazyWeakCompositions[n_Integer?NonNegative, 0, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[]] := LazyListCreate[{}, chunkSize, opts] lazyWeakCompositions[n_Integer?NonNegative, 1, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[]] := LazyListCreate[{{n}}, chunkSize, opts] lazyWeakCompositions[n_Integer?NonNegative, k_Integer /; k >= 2, chunkSize : _Integer?Positive : 10^5, opts : OptionsPattern[]] := With[{length = Binomial[n + k - 1, n]}, Module[{active = False, compositionsLeft = length, prevComposition = Join[ConstantArray[0, k - 2], {-1, n + 1}]}, LazyListCreate[ IteratorCreate[ ListIterator, (active = True) &, With[{realChunkSize = Min[chunkSize, compositionsLeft]}, If[realChunkSize === 0, {} ,(*else*) compositionsLeft -= realChunkSize; With[{taken = nextWeakCompositionsChunk[n, k, prevComposition, realChunkSize]},   prevComposition = Last[taken];   taken ] ] ] &, TrueQ[active] &, Remove[active, compositionsLeft, prevComposition] & ], chunkSize, opts, "Length" -> length, "FiniteList" -> True ] ] ] 
Fix typo in code.
Source Link
jkuczm
  • 15.2k
  • 2
  • 55
  • 89
Needs["CCompilerDriver`"] " #include \"WolframLibrary.h\" DLLEXPORT mint WolframLibrary_getVersion() { return WolframLibraryVersion; } DLLEXPORT int WolframLibrary_initialize(WolframLibraryData libData) { return LIBRARY_NO_ERROR; } DLLEXPORT void WolframLibrary_uninitialize(WolframLibraryData libData) {} DLLEXPORT int nextWeakCompositionsChunk( WolframLibraryData libData, mint Argc, MArgument *Args, MArgument Res ) { // number which compositions are calculated mint n = MArgument_getInteger(Args[0]); // length of composition mint k = MArgument_getInteger(Args[1]); // \" current \" composition tensor MTensor compositionT = MArgument_getMTensor(Args[2]); // actual data of \" current \" composition mint* composition = libData->MTensor_getIntegerData(compositionT); // one-dimentional, one-element tensor holding index of last non-zero // element in \"current\" composition MTensor lastNonZeroIndexT = MArgument_getMTensor(Args[3]); // pointer to index of last non-zero element in \"current\" composition mint* lastNonZeroIndex = libData->MTensor_getIntegerData(lastNonZeroIndexT); // number of compositions in returned chunk mint chunkSize = MArgument_getInteger(Args[4]); // dimensions of returned chunk tensor mint chunkDims[2] = {chunkSize, k}; // 2 dimentional tensor with chunk of compositions to be returned MTensor chunkT; libData->MTensor_new(MType_Integer, 2, chunkDims, &chunkT); // actual data of the chunk tensor mint* chunk = libData->MTensor_getIntegerData(chunkT); // index enumerating compositions in chunk mint i; // index enumerating elements of composition mint j; for (i = 0; i < chunkSize; i++) { // Moving to next composition does basically this: // {..., x, y} -> {..., x+1, y-1} if last element is non-zero, // {..., x, y, 0, ..., 0} -> {..., x+1, 0, 0, ..., y-1} // if last element is zero and y is last non-zero element. ++(composition[*lastNonZeroIndex - 1]); if (*lastNonZeroIndex == k -1) { --(composition[*lastNonZeroIndex]); } else { composition[k - 1] = composition[*lastNonZeroIndex] - 1; composition[*lastNonZeroIndex] = 0; } // If last element is zero it means we moved last non-zero element to // the left. If last element is non zero, then, obviously, it is the // last non-zero element. if (composition[k - 1] == 0) { --(*lastNonZeroIndex); } else { *lastNonZeroIndex = k - 1; } // Copy \"current\" composition to i-th row of chunk. for (j = 0; j < k; j++) { chunk[i*k + j] = composition[j]; } } // Return control over composition and lastNonZeroIndex tensors back to // Wolfram Language. libData->MTensor_disown(lastNonZeroIndexT); libData->MTensor_disown(compositionT); // Set chunk tensor as returned value of LibraryFunction. MArgument_setMTensor(Res, chunkT); return LIBRARY_NO_ERROR; } "; CreateLibrary[%, "weakCompositions"] ClearAll[nextWeakCompositionsChunk] nextWeakCompositionsChunk = LibraryFunctionLoad[%LibraryFunctionLoad[ %%, "nextWeakCompositionsChunk", {Integer, Integer, {Integer, 1, "Shared"}, {Integer, 1, "Shared"}, Integer}, {Integer, 2} ] 
Needs["CCompilerDriver`"] " #include \"WolframLibrary.h\" DLLEXPORT mint WolframLibrary_getVersion() { return WolframLibraryVersion; } DLLEXPORT int WolframLibrary_initialize(WolframLibraryData libData) { return LIBRARY_NO_ERROR; } DLLEXPORT void WolframLibrary_uninitialize(WolframLibraryData libData) {} DLLEXPORT int nextWeakCompositionsChunk( WolframLibraryData libData, mint Argc, MArgument *Args, MArgument Res ) { // number which compositions are calculated mint n = MArgument_getInteger(Args[0]); // length of composition mint k = MArgument_getInteger(Args[1]); // \" current \" composition tensor MTensor compositionT = MArgument_getMTensor(Args[2]); // actual data of \" current \" composition mint* composition = libData->MTensor_getIntegerData(compositionT); // one-dimentional, one-element tensor holding index of last non-zero // element in \"current\" composition MTensor lastNonZeroIndexT = MArgument_getMTensor(Args[3]); // pointer to index of last non-zero element in \"current\" composition mint* lastNonZeroIndex = libData->MTensor_getIntegerData(lastNonZeroIndexT); // number of compositions in returned chunk mint chunkSize = MArgument_getInteger(Args[4]); // dimensions of returned chunk tensor mint chunkDims[2] = {chunkSize, k}; // 2 dimentional tensor with chunk of compositions to be returned MTensor chunkT; libData->MTensor_new(MType_Integer, 2, chunkDims, &chunkT); // actual data of the chunk tensor mint* chunk = libData->MTensor_getIntegerData(chunkT); // index enumerating compositions in chunk mint i; // index enumerating elements of composition mint j; for (i = 0; i < chunkSize; i++) { // Moving to next composition does basically this: // {..., x, y} -> {..., x+1, y-1} if last element is non-zero, // {..., x, y, 0, ..., 0} -> {..., x+1, 0, 0, ..., y-1} // if last element is zero and y is last non-zero element. ++(composition[*lastNonZeroIndex - 1]); if (*lastNonZeroIndex == k -1) { --(composition[*lastNonZeroIndex]); } else { composition[k - 1] = composition[*lastNonZeroIndex] - 1; composition[*lastNonZeroIndex] = 0; } // If last element is zero it means we moved last non-zero element to // the left. If last element is non zero, then, obviously, it is the // last non-zero element. if (composition[k - 1] == 0) { --(*lastNonZeroIndex); } else { *lastNonZeroIndex = k - 1; } // Copy \"current\" composition to i-th row of chunk. for (j = 0; j < k; j++) { chunk[i*k + j] = composition[j]; } } // Return control over composition and lastNonZeroIndex tensors back to // Wolfram Language. libData->MTensor_disown(lastNonZeroIndexT); libData->MTensor_disown(compositionT); // Set chunk tensor as returned value of LibraryFunction. MArgument_setMTensor(Res, chunkT); return LIBRARY_NO_ERROR; } "; CreateLibrary[%, "weakCompositions"] ClearAll[nextWeakCompositionsChunk] nextWeakCompositionsChunk = LibraryFunctionLoad[%, "nextWeakCompositionsChunk", {Integer, Integer, {Integer, 1, "Shared"}, {Integer, 1, "Shared"}, Integer}, {Integer, 2} ] 
Needs["CCompilerDriver`"] " #include \"WolframLibrary.h\" DLLEXPORT mint WolframLibrary_getVersion() { return WolframLibraryVersion; } DLLEXPORT int WolframLibrary_initialize(WolframLibraryData libData) { return LIBRARY_NO_ERROR; } DLLEXPORT void WolframLibrary_uninitialize(WolframLibraryData libData) {} DLLEXPORT int nextWeakCompositionsChunk( WolframLibraryData libData, mint Argc, MArgument *Args, MArgument Res ) { // number which compositions are calculated mint n = MArgument_getInteger(Args[0]); // length of composition mint k = MArgument_getInteger(Args[1]); // \" current \" composition tensor MTensor compositionT = MArgument_getMTensor(Args[2]); // actual data of \" current \" composition mint* composition = libData->MTensor_getIntegerData(compositionT); // one-dimentional, one-element tensor holding index of last non-zero // element in \"current\" composition MTensor lastNonZeroIndexT = MArgument_getMTensor(Args[3]); // pointer to index of last non-zero element in \"current\" composition mint* lastNonZeroIndex = libData->MTensor_getIntegerData(lastNonZeroIndexT); // number of compositions in returned chunk mint chunkSize = MArgument_getInteger(Args[4]); // dimensions of returned chunk tensor mint chunkDims[2] = {chunkSize, k}; // 2 dimentional tensor with chunk of compositions to be returned MTensor chunkT; libData->MTensor_new(MType_Integer, 2, chunkDims, &chunkT); // actual data of the chunk tensor mint* chunk = libData->MTensor_getIntegerData(chunkT); // index enumerating compositions in chunk mint i; // index enumerating elements of composition mint j; for (i = 0; i < chunkSize; i++) { // Moving to next composition does basically this: // {..., x, y} -> {..., x+1, y-1} if last element is non-zero, // {..., x, y, 0, ..., 0} -> {..., x+1, 0, 0, ..., y-1} // if last element is zero and y is last non-zero element. ++(composition[*lastNonZeroIndex - 1]); if (*lastNonZeroIndex == k -1) { --(composition[*lastNonZeroIndex]); } else { composition[k - 1] = composition[*lastNonZeroIndex] - 1; composition[*lastNonZeroIndex] = 0; } // If last element is zero it means we moved last non-zero element to // the left. If last element is non zero, then, obviously, it is the // last non-zero element. if (composition[k - 1] == 0) { --(*lastNonZeroIndex); } else { *lastNonZeroIndex = k - 1; } // Copy \"current\" composition to i-th row of chunk. for (j = 0; j < k; j++) { chunk[i*k + j] = composition[j]; } } // Return control over composition and lastNonZeroIndex tensors back to // Wolfram Language. libData->MTensor_disown(lastNonZeroIndexT); libData->MTensor_disown(compositionT); // Set chunk tensor as returned value of LibraryFunction. MArgument_setMTensor(Res, chunkT); return LIBRARY_NO_ERROR; } "; CreateLibrary[%, "weakCompositions"] ClearAll[nextWeakCompositionsChunk] nextWeakCompositionsChunk = LibraryFunctionLoad[ %%, "nextWeakCompositionsChunk", {Integer, Integer, {Integer, 1, "Shared"}, {Integer, 1, "Shared"}, Integer}, {Integer, 2} ] 
Source Link
jkuczm
  • 15.2k
  • 2
  • 55
  • 89
Loading