2016-12-15 19:56:03 +0100 commented question Clean up recursive function Yes, that is the point I was trying to make. The function approximates the number of expected composites by dividing the range. But using floor implies the composite is at the end of each section and ceil implies it's at the beginning. Since they are both cumulative I'll probably use round, even though the others seem like a cleaner translation to the underlying density/distribution. Of course floor would give a good lower bound and ceil a good upper - as long as the error doesn't diverge. 2016-12-13 10:56:03 +0100 asked a question Clean up recursive function I've written a recursive function to approximate the count of k < K with 6*k +-1 both composite. It's accurate enough that now I want to: Clean it up* Use any theory/sage tools I should have started with Derive an equation, upper and lower bounds, etc. Write a paper on it I confess that I'm stumped by the recursion to equation/theory part. All my attempts at combinatoric, iterative, etc. failed. Can you help me clean it up and put it in it's best form? I may also have to do a print trace to help me see the math. I'm worried a code to math theory question won't fit anywhere. *Examples The line   ki = ceil(k/p)   If I use: ki = k/p     Treats k's as rationals not integers, range divisions always 'overlap'. ki = ceil(k/p)     k's are integers, but still divisions overlap too often. ki = floor(k/p)     Never overlap, but 'gaps' mean some k not counted. This is a coding detail that affects the underlying and any derived math. Function # twin composite counting function, recursive def twinc(k, flag, plist): # no more primes to try OR range < prime if len(plist) == 0: return 0 if k < plist[0]: return 0 p = plist[0]; del plist[0] ki = ceil(k/p) # half of twin is composite # ki = count of (k + pq): p | other half # try other primes where p misses if flag > 0: return ki + twinc(k - ki, flag, plist) # previous primes all missed, divide range into: # p | -1 or +1, p missed return twinc(ki, 1, copy(plist))*2 + twinc(k - ki*2, 0, plist) # end twinc  Test program pmax = 201 plist = prime_range(5, pmax) kmax = 5*max(plist) + 1 tlist = [(k, twinc(k, 0, copy(plist))) for k in range(5, kmax, 5)] P = list_plot(tlist, plotjoined=True) # direct plot fails, use stepping? P.show()  Update: Originally the flag was two bits for 0..3. While I was writing pseudo-code I realized that a boolean parameter can be eliminated with a pair of recursive functions. This looks closer to something that can be turned into an equation, yes? Alternate (pseudo-code) def twinCa(k, [p: plist]): ki = round(k/p) return k < p? 0: twinCa(k - ki*2, plist) + twinCb(ki, plist)*2 def twinCb(k, [p: plist]): ki = round(k/p) return k < p? 0: ki + twinCb(k - ki, plist)  2016-08-29 14:59:12 +0100 received badge ● Teacher (source) 2016-08-28 11:53:17 +0100 asked a question Finding lower bounds of a combinatoric function. I am working on a counting function involving prime numbers and managed to construct a sage function that is correct: def sum_X_subs(plist, plen): psum = 0 for k in range(plen-1): ptwo = 2^(plen-k) - 2 psum += sum(mul(psub)*ptwo for psub in Subsets(plist, k)) return psum # end sum_X_subs flist = [] plist = [5]; plen = 1 for p in prime_range(7, 41+1): plist.append(p) plen += 1 flist.append(sum_X_subs(plist, plen)) flist # [2, 52, 1162, 28196, 712250, 20379100, 696733730, 25016026280, 1042543611410, 47439135073960]  Now I need a lower bounds function for this: $$S = prime~range(5, n)$$ $$f(S) = \sum_{k=0}^{|S|-1} {\left(\left(2^{|S|-k}-2\right)\sum_{T\in \mathcal{P}_{k}(S)} {\prod_{t\in T} {t}}\right)}$$ None of the literature I've found covers anything close to it. (I have a vague memory of an Euler project problem that involved subset products). 2016-08-28 10:56:32 +0100 received badge ● Editor (source) 2016-08-28 10:53:26 +0100 answered a question Generating a LabelledBinaryTree from its string representation Ironically, the LabelledBinaryTree() command only lets you label the root. But you can apply a canonical labeling and then map those labels to the desired ones. import re words = re.compile('\w+') # keep int labels as ints def safeInt(word): try: return Integer(word) except TypeError: return word # labelled binary tree from string def lbt_from_str(lbstr): blank = words.sub('', lbstr) work = LabelledBinaryTree(blank).canonical_labelling() wstr = '%s'%work canon = [Integer(w) for w in words.findall(wstr)] labels = [safeInt(w) for w in words.findall(lbstr)] fixit = dict((canon[i], labels[i]) for i in range(len(canon))) return work.map_labels(lambda z: fixit[z]) # end lbt_from_str  This doesn't check if the string has a valid structure, and it expects all nodes to have labels. You could replace ' [' with 'None[' Testing: t = Permutation([1, 3, 2, 6, 5, 7, 4]).increasing_tree() print t show(t) s = '%s'%t test = lbt_from_str(s) print test show(test)  2016-05-24 21:50:02 +0100 answered a question Washed out colors when saving plots to file I think when the png is being constructed it doesn't handle the interaction between the opacity of the implicit plot and the point. Try adjusting your opacity. I also tested this: distance = 56; x,y,z = var('x y z') mypic = implicit_plot3d(x^2 + y^2 + z^2 == distance, (x, -10, 10), (y, -10, 10), (z, -10, 10), opacity = 0.5) for k in range(-9, 10, 2): mypic = mypic + point3d(vector([5, k, 4]), color = 'red', size = 10) mypic.save('mypic.png', compress = false) mypic  The inline image shows a line of six red points from up-left to down-right ( \ ), with a center gap where four are completely blocked by the sphere. The saved png shows all ten points in a line from down-left to up-right ( / ), with the center four slightly hazy. So the png is also swapping x with y for the points (the axis labels are correct). 2016-05-23 06:52:32 +0100 commented question How do I inform the sage webmaster regarding search results? Could someone post a wiki suggestion regarding redirects? Then my question can be de-wiki-fied. 2016-05-23 06:48:38 +0100 commented answer How do I inform the sage webmaster regarding search results? Thank you, I will make a habit of reporting both the old and new links to the site whenever I get redirected. 2016-05-23 06:48:31 +0100 received badge ● Scholar (source) 2016-05-23 06:48:26 +0100 received badge ● Supporter (source) 2016-05-23 02:08:42 +0100 received badge ● Student (source) 2016-05-22 04:50:33 +0100 received badge ● Organizer (source) 2016-05-22 04:16:37 +0100 answered a question Is there a functional way to manipulate cycles extracted from a digraph? With the help of a PARI programmer and OEIS regular I cleaned up my code: wief = DiGraph([prime_range(3600), lambda p, q: power_mod(p, q-1, q^2)==1]) sc = [[0]] + [sorted(c[1:], reverse=1) for c in wief.all_simple_cycles()] tbl = [sorted(filter(lambda c: len(c)==n, sc))[0] for n in range(1, len(sc[-1]))] for t in tbl: print 'n=%d:'%(len(t)), ', '.join("%s"%i for i in t)  and posted it to OEIS sequence A271100 Wieferich n-tuples. I used it to increase the number of listed terms from 15 to 300. We'll probably bump that to nearly 10k soon : ) *For the community, I'd like to add mouse-over hints describing what each piece of code does but I'm a noob. 2016-05-22 03:12:41 +0100 asked a question How do I inform the sage webmaster regarding search results? The Sage/Sagemath webmasters need to do the fix described below for all search engines: Some of the online sage resources that were moved still aren't visible to search engines yet. (Examples for Google MyWay, links trashed because I'm new) When I search for documentation: "sage intersection of lists" Most of the results are similar to: Sets — Sage Reference Manual v7.2: Basic Structures http://www.sagemath.org/doc/reference/sage/sets/set.html When I click the links I get: Redirecting ... http://doc.sagemath.org/html/en/reference/structure/sage/sets/set.html FIX: Google Search Help* says that the sage webmaster needs to: Create a "Google Webmasters Tools Account" - https://www.google.com/webmasters/tools/ then use the Fetch tool to get the spiders to you site so it can be indexed correctly - https://support.google.com/webmasters/answer/6066468?hl=en&ref_topic=6066464 2016-05-17 15:13:30 +0100 asked a question Is there a functional way to manipulate cycles extracted from a digraph? I wrote a worksheet that defines a digraph and extracts all of its simple cycles. pr = prime_range(3600) wief = DiGraph([pr, lambda p, q: power_mod(p, q - 1, q^2) == 1]) # could take a while! sc = wief.all_simple_cycles() sc.append([0, 0])  Then it uses mostly procedural code to: Remove the duplicate endpoints Reverse-sort the cycle Group cycles by length Sort the group of cycles Return the first (lexicographically earliest) cycle of each length maxn = 0 for c in sc: del c[0] # duplicate maxn = max(maxn, len(c)) c.sort(reverse=True) # reverse lexicographic # group by length tbl = [[c for c in sc if len(c) == n] for n in range(1,maxn)] n = 1 for t in tbl: # sort tuples by max element, get smallest tuple t.sort() print 'n=%d:'%(n), # why can't join handle a list of integers? print ', '.join(map(lambda m: '%s'%(m), t[0])) n = n + 1  The worksheet is adequate - it performs step 5 correctly, but a listing of the code will be posted as the PROG paragraph of an OEIS entry. So I would really like to reduce the numbered steps to as few lines as possible. Are there more graceful ways to tweak the output of one function before passing it on to another?